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ABSTRACT 


We investigate the behavior of acoustic cloaks and design an acoustic cloak applied to a thin 
spherical shell. We begin by examining the acoustic pressure field scattered from a spherical 
object, consider limiting impedance situations for a homogeneous sphere, and a thin shell. The 
second part of this thesis is focused on acoustic cloaking theory. Our research in this area con¬ 
centrated on the application of the transformational acoustics method as a solution to model and 
characterize physical parameters and fields involved in an acoustics cloak parameterization. A 
trade-off between an effective cloak and a possible realizable metamaterial dictates the charac¬ 
terization of the anisotropic inertia and stiffness parameters for the cloak design. Two limiting 
acoustic metafluids’ properties are explored. Lastly, we analyze the performance of the acoustic 
cloak when it is applied to scattering from a spherical shell. 
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CHAPTER 1: 
INTRODUCTION 


1.1 BACKGROUND 

Acoustic cloaking has been an area of intense research during the last decade. Being able to 
control the path of sound waves brings up a large number of possible applications, both in 
civil and military domains. New techniques in medicine, seismology and civil engineering are 
examples, as well as the new defense possibilities for submarines and surface warships. 

Nowadays, the search for a physical implementation of an acoustical cloak material, capable of 
redirecting incident acoustic radiation or mitigating scattered waves continues. Such a material 
is not a reality yet. Also, acoustical cloaking research is founded on acoustical scattering and 
echo signature prediction methods. There is not a unique solution that can be applied to every 
situation in which we may want to use acoustic cloaking materials. This limitation is intrinsi¬ 
cally related to the material properties and the shape of the target, the frequency and the sur¬ 
rounding fluid characteristics of the specific problem at hand. However, theoretical approaches 
have been providing very important conclusions and analytical, as well as qualitative, knowl¬ 
edge about the fundamental issues involved. Moreover, recent numerical work has provided 
notable insight about the limitations and potential for the development of acoustic cloaking for 
more interesting sorts of objects, materials and shapes. 

1.2 PURPOSE 

In this thesis, we begin our study by analyzing the acoustic field scattered from a homogenous 
sphere vibrating at a range of frequencies, in response to an incident plane wave. Initially, we 
consider the two limiting impedance boundary conditions: rigid and “pressure release.” This 
initial examination allows us to frame our research, as well as gain a better understanding of 
spherical scattering phenomena. Subsequently, we will extend the scattering analysis to a fluid 
loaded thin spherical shell and again consider scattering for a range of frequencies that includes 
the shell resonances. 

Subsequently, we develop an acoustic cloak analytical solution applied to the thin spherical 
shell. Finally, we proceed with several computational simulations, varying the incident fre¬ 
quencies and the characteristic parameters of the thin shell. We perform these simulations with 
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and without the application of the acoustic cloak in order to evaluate the scattering resonance 
frequencies, amplitudes and the back scatter reduction we are able to achieve when we apply an 
acoustic cloak to the system. 

1.3 RELATED RESEARCH 

Acoustics cloaking became a hot topic of research over the last few years. Former work began 
to exploit the theory of invisibility applied to electromagnetic waves. The research on imple¬ 
mentation of a cloaking material has been spurred by work in transformation of coordinates, 
in both electromagnetic and acoustics fields. In 1996, Pendry el al. [1] developed a coordinate 
transformation system that would allow for adaptive scaling of the electromagnetic parameters 
of permeability and permittivity as necessary, while performing studies of Maxwell’s equations 
in complex geometries. In 2003, Greenleaf et al. [2] present a notable result in studies on con¬ 
ductivity that indicate the possibility of cloaking. They developed anisotropic conductivities 
in 3D that could provide analogous measurement results in voltage and current measured by a 
body with homogeneous and isotropic conductivity. While most of research has been focused on 
developing electromagnetic cloaks, Leonhardt [3] postulated that cloaking theory might apply 
equally well to acoustic waves. In 2007, Cummer and Schrurig [4] provided a 2D transforma¬ 
tion type for acoustics waves, using anisotropic mass density and showed that, for this particular 
case, it is theoretically possible to achieve an acoustic cloak at fixed frequencies. This contra¬ 
dicted a contemporary result, presented by Milton et al. [5] where it was argued that the acoustic 
equations of motion are not form invariant, and consequently it would not be possible to employ 
an analogous acoustic cloak with the theory developed by Pendry et al. [1] in electromagnet¬ 
ics. Chen and Chan [6] pursued the line of research based upon a mapping of the 3D acoustic 
equations similar to those used in conductivity studies, introducing a spherical-Bessel functions 
expansion method. Cummer et al. [7] stated the necessary constraints for anisotropic mass den¬ 
sity and radial bulk modulus dependence in their derivation, and hypothesized new methods of 
implementation for an acoustics cloak with no direct relation or comparison with electromag¬ 
netics. The present work follows more closely the work developed by Norris [8], although with 
different conceptualization of the original and transformed regions. Norris developed a trans¬ 
formational acoustics method using a change of variables and a one-to-one mapping between 
two coordinate references. The one-to-one mapping is applied throughout the cloaking region, 
with the exception of a single point or surface that constitutes the inner boundary of the cloaking 
region. The definition of an Inertial Cloak (IC) and the metafluids necessary for its construction 
is provided. The IC is characterized by anisotropic density, and Norris exposes the problem of 
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an infinite mass density for a perfect cloaking using IC metafluids and provides the solution: a 
combination of IC metafluids and materials that are anisotropic in their stiffness, referred to as 
Pentamode materials (PM), already proposed much earlier in the literature [9]. This solution 
predicts the possibility of achieving an acoustic cloak material with realistic physical properties. 

1.4 THESIS OUTLINE 

Chapter 2 addresses the scattering phenomena observable for a homogeneous sphere and for a 
thin spherical shell immersed in a fluid. Different qualitative characteristic impedance values 
are considered, as well as the boundary conditions for each case. The total and the scattered 
pressure fields are analytically determined, and the radial displacement and particle velocity 
variations are investigated. Additionally, the scattered energy is evaluated for the prescribed 
cases of the acoustically hard/soft spheres and the thin spherical shell. 

Chapter 3 is dedicated to the development and application of the transformational acoustics 
technique. Coordinate transformations are examined and the different possible solutions, as 
well as the inherent trade-offs, are described. We apply the selected analytical solution to the 
spherical thin shell and analytically evaluate the scattered pressure amplitude and energy when 
an acoustic cloak is applied to the shell. 

Chapter 4 presents computational results performed for each analytical model and extends the 
results achieved in the previous chapters, allowing us to extract more complete solutions and 
conclusions about the theoretical solutions. 

Chapter 5 briefly summarizes the main points of the work and enumerates the main conclusions. 
Additionally, further research topics and considerations are suggested for continuation of this 
research effort. 
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CHAPTER 2: 

ACOUSTICS SCATTERING 


This chapter presents the fundamental concepts related to the scattering problem. We begin 
our study of the scatted pressure field by performing the mathematical development for the two 
limiting cases of the acoustical impedance for a sphere: rigid and “pressure release” boundary 
conditions. In what follows, the analytical scattering approach is extended to a thin spherical 
shell object. 

2.1 PROBLEM DESCRIPTION 

We consider a submerged sphere and an incident plane wave traveling along the z axis, in the 
positive direction, as illustrated in Figure 2.1. 



Figure 2.1: Incident plane wave on a spherical obstacle 


A spherical boundary constitutes a very special shape that can be examined analytically. The 
reason is that it is possible to easily define this boundary in an orthogonal coordinate system 
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(spherical coordinates, in this case), by setting one of the coordinates to a given value. More¬ 
over, the wave equation is separable in this coordinate system [10]. The pressure field turns out 
to be represented by a linear combination of radial modal amplitudes and spherical wave har¬ 
monics. Furthermore, the spherical axisymmetry leads to the non-dependence of the azimuthal 
acceleration distribution. That is, there is no azimuthal angular dependence (f), and consequently 
the spherical wave harmonics are reduced to Legendre polynomial expressions. 

The scattered pressure is defined as the change in a sound field due to a boundary being intro¬ 
duced into the fluid [10]. For our purposes, we consider that this “introduced boundary” to be 
the spherical shell to which we will apply an acoustic cloak. 

The total pressure field pr, from the outer surface of the shell’s point of view, is considered to 
be the sum of an incident and scattered pressure field: 

Pt = Pi+Ps (2.1) 


2.2 INCIDENT AND SCATTERED PRESSURE FIELDS 
FROM A SPHERE 

We start our study by defining the incident and the scattered pressure fields for a homogeneous 
sphere, considering the 3D Helmholtz equation in spherical coordinates and finding its solu¬ 
tions: 

(V 2 + k 2 )p = 0 (2.2) 


where 

v 2 „= . 1 9 ( Aae a ’’) + 1 d '‘r 

r 2 dr V dr ) r 2 sin 6 dd V 89/ r 2 sin 2 6 d 2 (j) 


(2.3) 


is the Laplacian operator in spherical coordinates. The pressure field is denoted by p, k is the 
wave number k = u/c (m -1 ), c o is the angular velocity (rad/s) and c is the speed of sound 

(m/s). 


The solution of the Helmholtz equation is determined by applying the method of separation of 
variables, and results in an expression for the incident and the scattered pressure fields in terms 
of products of spherical Bessel functions and associated Legendre polynomials. 


2.2.1 INCIDENT PRESSURE FIELD 

The incident pressure field is taken to be a plane wave traveling along the z axis as illustrated in 
Figure 2.1: 

p r = Pie ikz e~ iwt (2.4) 
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where e lwt denotes the traveling wave’s time dependence and P, is the incident pressure am¬ 
plitude. 

We expand this expression in terms of a spherical basis of solutions of Equation (2.2) [11]: 

Pl = p ie ^r cos e) e -iwt 

= i c mnjn(kr) + d mn n n (kr)\P™( cos ey m *e- iwt (2.5) 

n m 

where j n (kr) and n n (kr) are the spherical Bessel’s functions of first and second kinds, respec¬ 
tively, and c mn and d mn are modal scaling coefficients of the pressure amplitude, associated with 
their corresponding spherical Bessel functions. P m ' p denotes the azimuthal angular dependence. 

The incident pressure is, by definition, bounded (finite) at the origin of the coordinate system, 
which we take to be the target center of the sphere. Therefore, the term given by the spherical 
Bessel’s equation of the second kind, n n (kr) (which —>■ oo as r —>• 0), is unphysical and we must 
set d mn = 0. The coefficient c mn = (2 n + 1)(/’)'' is derived simultaneously with the derivation 
of the definition of the plane wave (Equation (2.4)) in spherical coordinates (Equation (2.5)). 

The incident pressure field is then given by: 


Pi(r, 0) = Pi ^(2 n + l)i n j n (kr )P n (cos9)e lwt (2.6) 


2.2.2 SCATTERED PRESSURE FIELD 

As the incident pressure field satisfies the Helmholtz equation, so too does the scattered pressure 
field ps- 

A general solution includes incoming and outgoing traveling pressure waves. Consequently, 
spherical Bessel functions of the third kind are used with corresponding coefficients: 

PsM.f) = P.J2Y1 i e ™ h ™( kr ) + / m „hW(^)]^ m (c 0 s9)e ,m *e-" (2.7) 


where /( is scaled, in this equation, by e mn and The functions hi!:., and h.!.: are the 

spherical Bessel’s functions of the third kind, and e mn and f mn are, similarly to c mn and d mn , 
the modal scaling coefficients corresponding to their corresponding spherical Bessel functions. 

The spherical Bessel functions of the third kind are also called spherical Hankel functions, and 
are appropriate for the description of spherically traveling waves. Having defined an e~ lwt time 
dependence, hmn and hmn denote outgoing and incoming spherical waves, respectively. 

We include only the spherical Hankel functions of the first type, hmn, since we are interested in 
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outgoing traveling spherical waves, and the expression forps, becomes: 


Ps(r,9 ) = P s [e mn h (p n (kr)]P n (cos(9))e lwt (2.8) 


The scaling coefficient e mn is determined once the rigid or the “pressure release” boundary 
conditions at the surface of the sphere are imposed. 


2.2.3 Rigid Boundary Condition 

The rigid boundary condition is enforced by requiring that the normal component of the particle 
velocity vanishes at the surface of the sphere, r = a [10], and require: 


dp T 

dr 



0 


(2.9) 


We shall denote, by considering Equation (2.1), that the imposition of the rigid boundary con¬ 
dition requires: 


dpi dp s 
dr dr 


( 2 . 10 ) 


Because the incident plane wave is independent of 0 and the obstacle has spherical symmetry, 
the coefficient e mn , in Equation (2.8), is found to be independent of m, and is found to be: 

e m „ = e n = —(2 n + 1 )i n J t ^ a \ (2.11) 

hp ( ka) 


The pressure field scattered from a rigid boundary (pf) is, therefore, given by: 


PsM) = [( 2W + l) r s tl n ) ( kr ) Pn{cOs(9))e 

n tin \rCCL) 


—iwt 


( 2 . 12 ) 


For k/r —» oo the outgoing spherical Hankel function has the asymptotic form [10] [12]: 

h£\kr) « (2.13) 

kr 

and n corresponds to the (finite) number of the normal mode. Substitution of Equation (2.13) 
in Equation (2.12) yields the far-field scattered pressure as a function of 6 and r: 



i Pi c i(kr) 
kr 


Y, (2n + 1) 

71=0 


jn(ka) 

h^\ka) 


P n (cos(9))e~ iwt 


(2.14) 
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Equation (2.14) depicts the far-held pressure due to scattering of a plane wave from a rigid 
spherical obstacle. 

2.2.4 “Pressure Release” Boundary Condition 

An ideal “pressure release” boundary condition requires that the total pressure, given by Equa¬ 
tion (2.1) vanishes at the surface of the sphere [10]. That is: 

Pi+Ps = 0 (2.15) 


at r = a. 

The coefficient e mn , in Equation (2.8), is derived using the same approach as for the rigid 
boundary case. Likewise, this coefficient does not depend upon m: 

e mn = e n = —(2 n + 1 )i n (2.16) 

hn (ka ) 


The scattering pressure held is, therefore, expressed as: 

■n 3n(ka) 


Ps (r, 0) = -Pi ( 2n + l ) 1 


hn\ka) 


n =0 


h£\kr) P n (cos9)e 


—iwt 


(2.17) 


The far-held condition stated in Equation (2.13) is, again, considered in order to obtain the the 
far held “pressure release” scattered pressure held: 


Ps„(.r, 6) = f; (2n + \)^PPp n (cosO)e.- M 


( kr ) 


hn\ka ) 


n=0 


(2.18) 


2.3 SCATTERED PRESSURE FIELD FROM A THIN 
SPHERICAL SHELL 

Equations (2.12) and (2.17) describe the scattered pressure helds, associated with the two lim¬ 
iting mechanical impedance conditions that lead to the acoustically rigid and “pressure release” 
boundary conditions. 

We now determine the pressure held scattered from a submerged, evacuated, thin spherical 
shell. 

The shell is considered to be a thin shell if the thickness (h) is less than A, s /20, where A s is the 
shear wavelength, and is determined by considering the shear velocity, c s = \J G/p s , where 
G = E/ 2(1 + v) is the shear modulus, E is the Young’s modulus, v is the Poisson’s ratio, and 
p s is the density of the structural material. 
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The incident pressure field is considered to be a plane wave traveling in the positive z direction, 
as in the previous discussion. The former analysis lets us anticipate that the scattered pressure 
field will be expressed in terms of spherical Bessel functions of the third kind and spherical 
harmonics. However, in this case, the scattered pressure is not defined by a limit boundary con¬ 
dition, and will require a slightly different set of considerations. Nonetheless, the total pressure 
is still given by Equation (2.1). Substituting Equations (2.6) and (2.8) in to this equation, we 
obtain the general mathematical expression for the thin spherical shell case: 

OO 

Pt S = ^ J}( 2n + l )PUkr) + R n h£\kr)]P n (cQS0)e- iwt (2.19) 

n=0 


and 

OO 

P? S (a, 9) = P i J2Pn S Pn(cos9)e~ iwt (2.20) 

n =0 

where p^ s represents the nodal amplitude of the total pressure amplitude at r = a. 

The scattered pressure field from a thin spherical shell is, consequently, defined by the expres¬ 
sion: 

OO 

P T s S = p i J2l R nh { n ) (kr)\P n (cos9)e- iwt ( 2 . 21 ) 

n =0 

where R n is the modal coefficient associated with the Hankel function, hn \ and determines 
the amplitude of the outgoing scattered waves, similar to the c mn and e mn coefficients in the 
previous section. 

We now need to apply the equations that determine the R n coefficients in order to define the 
final mathematical expression for p^ s . We begin by calculating the amplitude of the radial 
displacements and stresses of the fluid at the surface of the spherical shell corresponding to the 
incident pressure field. 

We assume that the normal velocity and displacement of the system fluid/submerged spheri¬ 
cal shell match, in amplitude and phase, at the outer surface of the shell. That is, there is no 
cavitation between the fluid and the shell. This assumption presumes small amplitude vibra¬ 
tions of the system, which is an adequate assumption given the difference in the characteristic 
impedances of the two media. 

In order to determine the displacement of the outer surface of the shell we exploit the equa¬ 
tions of motion for a spherical shell as they are defined in [10], where it is assumed that the 
spherical shell is thin enough and the frequency is low enough for it to be possible to “ignore 
flexural stresses as compared with membrane stresses.” Additionally, no torsional axisymmetric 
motions are taken into consideration, and so the tangential component of motion in 0 vanishes. 
The functions w and u represent the mid-surface radial and tangential motions of the shell, re¬ 
spectively, and are independent of 0. The elasticity equations for w and u, approximated for the 
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thin material shell, are: 


0 2 u 


Ou 


(i + 0 ) taa + cot 9— - (is + cot 0)u 


109 2 

n 0 2 VJ 

-P + 


09 


-0 


,0 3 w 


(1 + v) + f3 2 (u + cot 2 9) 


09 3 
Ow a 2 u 

- -7T- =0 


J 09 


( 2 . 22 ) 


and 


^ + 2f m te 92u 


09 3 

+ cot 9 


09 2 


(1 + i/)(l + (3 2 ) + /3 2 cot 2 0 


Ow 

00 


(2 - v + cot 2 0)/3 2 - (1 + v) 


n n O^VJ 2 (9 3 w 

u - 0 - 2/3 cot 0- 


<90 4 


<90 3 


+ /3 2 (1 + z/ + cot 2 0)^-^ - /3 2 cot 9(2 - i/ + cot 2 0)^ 
<70 z at/ 

, . a 2 w) (1 — z/ 2 )a 2 

- 2(1 + - -s- = -^—Eh 


(2.23) 


where a is the neutral-surface radius of the shell, h is the thickness of the shell, (3 = \J 00 and 
v and E are the Poisson’s ratio and the Young’s modulus, respectively, [13]. The parameter c p 
is defined as the low-frequency phase velocity of compressional waves in an elastic plate [10]: 


* s <2 - 24) 
and finally, the total pressure pr is given by Equation (2.1). 

The radial and tangential displacements u and w can be derived in terms of the Legendre poly¬ 
nomials of the first kind [10]: 


where r/ = cos 9. 


u (v) = 

n =0 


w{rj) = 


TT A - 2 dP n 

U W1 v drj 

(2.25) 

oo 

Y. W»P n (ri) 

n=0 

(2.26) 


Substitution of Equations (2.25), (2.26) and (2.20) for the total pressure evaluated at r = a 
in Equations (2.22) and (2.23), allows us to conclude that the coefficients U n and W n must 
satisfy [10]: 


[Q 2 - (1 + p 2 ){u + X n - 1 )}U n - [(3 2 {u + A n — !) + (! + v)]w n = 0 (2.27) 
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and 


— X n [f3 2 {y + A n — 1) + (1 + v)\U n 

+ [ft 2 - 2(1 + is) - (3 2 X n (is + X n - 1 )]W n = - a {1 ~ h U ) p^ 5 (2.28) 


where ft 2 = aco/c p and A n = n(n + 1). The remaining variables have the same definition as in 
the previous equations. 


We use Equations (2.27) and (2.28) to obtain W n (U n can be obtained in the same way): 

a 2 (l - is 2 ) 

r iPn 

ft 2 — ya 


W — — " _ ~—Lp. v TS 

Eh lVn 


(2.29) 


12 4 — 12 2 (7 a + 2e + /3' 2 X n a) — X n ((3 2 a + e) 2 + 2 ^: 70 ; + /3 2 a 2 A n 7 
with a = is + X n — 1, 7 = 1 +/5 2 and e = 1 + is. 

^ A A 

Let d = wf + u6 + i'(j) denote the the vector displacement of a fluid element with density p/ at 
the surface of the spherical shell. Then the linearized Euler’s equation [13]: 


d 2 d 

pf W ~ 


-Vp ? 5 


(2.30) 


can be used to relate W n to the total modal pressure. The parameter Vpr is the gradient of the 
total pressure at the surface of the sphere. 


We reduce Equation (2.30) to its respective radial components, to find the modal amplitudes of 
the total pressure gradient: 


d 2 w n _ d( pt S ) u 

dt 2 dr 


(2.31) 


The modal specfic acoustic impedance z n of the spherical shell is defined as: 


= P f 


P 


TS 


W„ 


(2.32) 


Given the time harmonic dependence W n = —iuW n , and substituting Equation (2.29) in Equa¬ 
tion (2.32), we obtain an expression for the modal specific acoustic impedance for the spherical 
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shell: 


Z n 


iPsCp h 
f2 a 

f2 4 — Vt 2 (^a + 2e + /3 2 X n a ) — A n (/3 2 a + s) 2 + 2^70: + (3 2 a 2 A„7 

f2 2 — 7 a 


(2.33) 


The coefficient R n is determined by calculating the derivative of the total pressure (Equation 
(2.19)): 

o j'S 00 

= y^[(2n + l)i"7cjp(/cr) + J R n /i^ (1) (/cr)]P n (cos6>) (2.34) 

n =0 

Substituting the modal coefficient evaluated at r = a into Equation (2.31) and eliminating W n 
from Equation (2.32) we obtain the scattering coefficient R n at the outer surface of the thin 
spherical shell (r = a): 


R 


n 


— (2n + l)i n 


Znj'n(ka) 

z n hn l \ka) 


izgj n (ka ) 

iz a hn\ka) 


(2.35) 


where the z a = pfCf (the specific acoustic impedance of the fluid) has been introduced. 


Substituting Equation (2.35) into Equation (2.19) we obtain the complete expression of the total 
pressure field due to a plane wave scattering from a thin spherical shell: 


ts \ Vo 1 11 ' n T • 11 \ z n kj (ka) iZajnika) 

Pt =2^( 2n + 1 ) z Jnikr)- , (1) x - . 77TT77 h n\kr) 


n=0 - z n khn L> {ka) - iz a hn\ka) 

and, using Equation (2.21), the scattered pressure is 


P n { cos 9)e 


—iwt 


(2.36) 


pf = -Pi £>n + l)t" - h^(kr)P n (cosO)e~ iu!t (2.37) 

z n hn ( ka ) — iz a hn {ka) 


n=0 
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CHAPTER 3: 

TRANSFORMATIONAL ACOUSTICS 


This chapter addresses the transformational acoustics method as an approach to design an acous¬ 
tic cloak. We consider an original coordinate system, where the spherical shell, of radius a is 
immersed in a cloak material, designated in general as a metafluid. The anisotropic properties 
of the acoustic cloak lead to a definition of the pressure field within the cloak as a ’’pseudo¬ 
pressure field” [8], proportional to the average of the compressive stress but not defined by the 
regular wave equation. The transformed coordinate system makes the spherical shell appear 
smaller, as the transformed radius has length 5 < a. Figure 3.1 illustrates the original and the 
transformed coordinate systems. We can observe the decrease in radius length as a net effect 
of the transformational acoustics application and that the exterior surface of the acoustic cloak 
remains the same (r — b). In the transformed coordinate systems, the metafluid is governed by 
the regular wave equation and is, therefore, homogeneous and isotropic. 



Figure 3.1: Transformational acoustics: (a) original and (b) transformed coordinate systems 


An acoustic metafluid is a material characterized by peculiar physical properties that do not exist 
in nature. Therefore, they would have to be manufactured to possess these unusual features. In 
theory, such material cloaks simulate a coordinate transformation where the wave equation of 
the pressure field within the acoustic cloak mimics the regular wave equation of the pressure 
field outside the cloak. This result is the scattering obstacle appearing to have shrunk to a 
smaller size. Indeed, ideally it would be zero and the obstacle would be shrunk to a single 
point which does not scatter energy, corresponding to a perfect cloak. However, as it will be 
shown, this ideal limit is not physically possible to achieve. Consequently, there will always 
be a scattered pressure wave. Nonetheless, scattering can be mitigated and a portion of the 
incident wave is bent around the object to be cloaked as a net effect of the physical properties 
of the metafluid. 

The most general of the metafluids contains anisotropic properties in mass density and stiff¬ 
ness. These metafluids are called Pentamode-Inertial Cloak (PM-IC) metafluids. We begin by 
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exploring the fundamental physical principles and the corresponding mathematical implemen¬ 
tation of a pure Inertial Cloak (IC), which is based on the imposition of anisotropy in the mass 
density of the fluid alone. Secondly, we explore the characteristics of the Pentamode materials 
(PM), considering the principles of elastodynamics applied to isotropic solids, ideal fluids, and 
extending these concepts to a pure PM material which is characterized by isotropic mass den¬ 
sity and anisotropic stiffness. The general equations applicable to PM-IC metafluids are also 
developed in order to describe the more general problem and allow for conclusions about the 
strengths and limitations of the cloak implementations. Lastly, we examine the application of 
pure IC and pure PM cloaks to a spherical geometry. 

We will consider an ideal metafluid surrounding a spherical shell where the anisotropic proper¬ 
ties are manifest and vary continuously. 


3.1 ACOUSTICS CLOAKING THEORY 

The conditions for the “invisibility” of the cloak require that the acoustic pressure field and the 
particle velocity of the external fluid and the acoustic cloak must match across the outer surface 
of the cloak. In other words, there is no change in acoustic impedance between the fluid and 
the outer surface of the cloak. Within the cloak, the acoustic pressure field is called “pseudo¬ 
pressure” and mimics the acoustic pressure waves of the isotropic and homogeneous exterior 
fluid [8]. 


The acoustic transformation maps the “pseudo-pressure” wave equation into the regular wave 
equation, given by: 


V 2 X P- 


1 d 2 p 
c 2 dt 2 


0 


or 


KO _ 2 


= 0 

Po dt 2 


(3.1) 


where Vx denotes the differentiation with respect to the coordinates {A",, X 2 , X :i }. p 0 is the 
mass density, /Co is the bulk modulus of a homogeneous and isotropic fluid, and we note that 
the phase speed c can also be defined by c 2 = /C 0 /po for adiabatic acoustic processes of small 
amplitude, in fluids. 


The general (for a PM-IC) “pseudo-pressure” wave equation will be derived in this chapter. It 
is given by: 


/CQ : V a 


Q 'QVjc/J 


d 2 p 
dt 2 


= 0 


(3.2) 


where p is the pseudo-pressure within the metafluid, and /C has the dimensions of stiffness, 
which will be defined specifically for each fluid. The tensor Q will be defined in detail in the 
PM materials section, and the operator (:) can be defined in terms of the cartesian components 
of the tensors A and B by A : B = Cfp = AijBij, where in the last step we have used 

* j 

the Einstein summation convention in which repeated indices are summed. 


16 



The transformational acoustics method that we will follow defines Q as a divergence free but 
non-unique stress tensor [8]. Indeed, fundamental characteristics of the tensor Q are a departure 
from the theory presented in [5], where an isomorphic relation between all the equations and 
field variables is required. That is, all fields and variables are mapped, in a one-to-one basis, 
from an original region to a transformed region. Such an approach would lead to a harmonic 
coordinate deformation which is ineffective for the purposes of cloaking [8] [5] [14]. 


Under the previous assumptions, the equivalent linearized displacement equations of motion are 
imposed: 

d 2 d 

^ _1 V a ,p+ —= 0 (3.3) 


and 


<7 = Q 


dt 2 


(3.4) 


where g is the second order tensor mass density of the fluid, d is defined as in Equation (2.30) 
and a defines the stress tensor of the metafluid. These two conditions are employed in the 
following sections, in order to define the dynamic equations for IC and PM acoustic cloaks. 


We now define notation for the transformation of coordinates and the definition of the variables, 
vectors and tensors elements. In general, and if not mentioned otherwise, we shall use lower 
case letters to refer to the original orthogonal coordinate system and capital letters to correspond 
to a transformed orthogonal coordinate system. 

We let x refer to the original Cartesian coordinate system {xi, x 2 , x 3 }, where we assume a 
cloaked region {r : a < \r \ < b}, where r = a is inner surface of the cloak and r = b is the outer 
surface of the cloak. The transformed coordinate system X with components {Xi,X 2 ,X 3 }, 
will allow us to obtain a description of the “pseudo-pressure” field satisfied in the x coordinate 
by a standard wave equation in the X coordinates. 

Additionally, in the component form decomposition, V 3 , • A = if A is a vector, and 
V* • A = C -^ L if A is a tensor. 

We define unit vectors x = x/r and X = X/R in the original and new coordinate systems, 
correspondingly, and the transformation is performed by a one-to-one mapping. 


The deformation tensor F is given by: 


F = W x x (3.5) 

with inverse: 

F 1 = V x X (3.6) 

J is the Jacobian of the coordinate transformation, and is defined as J = det (F) = |F|. 
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The deformation tensor F can be written, alternatively, using a polar decomposition: F = 
V/ Ra , where Rl is an orthogonal matrix and V l is the symmetric positive definite left stretch 
tensor that satisfies FF 7 = V|. 

In practice, the V L and H L tensors can be seen as compression and rotational tensors within 
transformation acoustics, respectively. 

The variation of a point element and a volume element, between coordinate systems is given 
by: 


dx = Fdx 


(3.7) 


and 


dv = \dx 3 • (dxi x dx 2 )\ 

= |F|(9X 3 • dX x x dX 2 ) 

= \F\dV 

dv = JdV (3.8) 


Denoting ds, dS, as the element areas and n and TV. the correspondent unit normal vectors, of 
both coordinates systems, we have from Nanson’s formula, [15] [8]: 

h • (j _1 F)ds = NdS (3.9) 


Employing Gauss’s Theorem [16]: 



Wx-FOV = 


Bo 


F TJ N,dS 


9B 0 


(3.10) 


♦ 

where B 0 and 9B 0 denote a volume and a surface area in the X coordinate system. Letting, 
F u = <5/j, we have V^-F = 0 and F uN I = Nj 


Combining results obtained by Equations (3.9) and (3.10), we derive: 


0 = 



n*(J _ F)ds 


dB 


(3.11) 


Employing again the Gauss Theorem: 


0 = / / n*(J _1 F)ds = / / / V x -(J~ l F)dv 
J J dB J J JB 


(3.12) 
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where B and OB denote an arbitrary volume and its bounding surface in the x coordinate system. 
Therefore: 

W x -{J~ 1 F)=0 (3.13) 

The above result will be used in the coordinate transformation from x to X as well as provide a 
candidate tensor for designing the pentamode stiffness matrix. 


3.2 INERTIAL CLOAK 

Using the coordinate transformation (Equation (3.5)), the polar decomposition of the deforma¬ 
tion tensor F, and the relationship Equation (3.13), Norris [8] proves the relationship: 

V 2 xV = JV X • ( J-W 2 X x p) (3.14) 

Specifically, the terms involving the acoustic pressure within the metafluid on the right hand 

♦ 

side of Equation (3.14) are equivalent to the Laplacian of the pressure in the X co-ordinate 
system. 

In order to determine the IC density and stiffness parameters, we consider the fundamental 

—► 

relation between the pressure and fluid displacement vector d [10]: 

p = -JCW x -d (3.15) 

as well as Euler’s linearized Equation (3.3), rewritten here for convenience [10]: 

pB d 

Q~ 1 ^xP+-gp= 0 (3.16) 

The pure IC fluid is considered to be inviscid, with isotropic bulk modulus /C and anisotropic 
density q (represented by a symmetric tensor of second order). Taking two time derivatives of 
Equation (3.15) and the gradient of Equation (3.16) and combining results to eliminate d, we 
obtain: 

K,,V X ■ {(f { V xP ) = (3.17) 

For the left hand side of Equation (3.17) to match the right hand side of Equation (3.14) we 
assign the stiffness and mass density parameters of an IC to be respectively: 

K. = JIC 0 (3.18) 

0 = p o J(X~ 1 ) 2 (3.19) 

where /Co and p 0 are the outer fluid constants bulk modulus and density, respectively, and J and 
V are related to our coordinate transformation. 
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3.2.1 Boundary Conditions for an Inertial Cloak 

As mentioned in the previous section, the pressure field and the normal velocity must match at 
the outer surface of the cloak. Explicitly, continuity of the normal stress (or pressure) p requires 
the same amplitude and phase values whether we are using Equation (3.1) at the exterior fluid 
side r = b + , or Equation (3.17) on the IC side r = b~, that is, \p] = p{b + ) — p(b~ ) = 0 at r = b. 

The outer surface of the cloak and the exterior fluid surface area contiguous with the cloak must 
be identical 1 : ds = dS 1 . 

We employ Equation (3.9), rewritten here for convenience, and Equation (3.19), at r = b: 

h • (j“ 1 F)ds = NdS 
(j -1 F T )n = N 
J- 1 FF r n = F N 

J- l \ 2 n = F N (3.20) 

and (po£ _1 ) = J _1 V 2 . The last line of Equation (3.20) is a purely kinematic condition [8], and 
ensures conservation of mass at r = b. 

Additionally, we consider Equation (3.16) and note that: 

d 2 d / \ 

^'^2 = (3-21) 

Similarly, exterior to the cloak: 

d 

= ' (do^VccP) (3.22) 

Therefore, combining Equations (3.21) and (3.22) we obtain the relation: 

h-g^Vp = pa~ l h-V x p (3.23) 

which is the requirement for continuity of the normal displacement at the outer boundary of the 
cloak. 

So far, we have discussed the governing equations of the IC metafluid, and the boundary con¬ 
ditions at the outer surface of the cloak. We stated previously that the IC is characterized by 
an anisotropic mass density, which manifests itself by our requirement p = p 0 J(V -1 ) 2 , and 
where (V -1 ) 2 is not proportional to I. The coordinate transformation that allows us to model 
the metafluid IC is equivalent to a contraction without rotational deformation. The mapping be¬ 
tween the original coordinate system and the deformed coordinate system can become singular 

'Small amplitude vibrations are assumed, similarly with conditions assumed in Chapter 2. 
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if the surface of the scatterer is transformed to a single point. In this instance, the density of the 
metafluid becomes more extreme, and finally goes to infinity when the surface is transformed 
to a point. [8]. 


3.3 PENTAMODE MATERIALS 

The PM material constitutes an alternative to the IC metafluid. 

A PM material is so-named because it has only one independent elastic modulus. Indeed, 
pentamodes are a generalization of isotropic fluids for which the only elastic eigenmode is 
hydrostatic stress, or pure pressure. 

We define the governing equations for a PM and verify that the “pseudo-pressure” field also 
satisfies the regular wave Equation (3.1) upon transformation as was done for the IC. We begin 
the implementation of the fundamental equations for a PM by considering the 3D version of 
Hooke’s law [8], given in the form: 

&ij C'ijkl^-kl (3.24) 

in Cartesian coordinates, where a is the stress tensor, C is the fourth order stiffness tensor and 
e is the strain tensor, defined as [17]: 


e 


ij 


1 r dui 

2 .dxj 


+ 


duf 

dxi- 


(3.25) 


Because they are symmetric, stress and strain tensors can be conveniently written as [6 x 1] 
order vectors, and the stiffness tensor C can be expressed as a square [6 x 6] matrix 2 : 



/ 


( en ') 


022 


^22 

J 

033 

} 

e 33 

a = 

V / 20‘23 

> e = 



V2<7 31 


V^^31 


V V2a u j 


V V2e 12 / 


2 Matrix notation used in this thesis follows notation used in [8], which diverges from the traditional notation 
[18] [19]. The notation we are using has the appeal that matrices are truly vectors and tensors, whereas they are 
not in the conventional matrix notation. 
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and 


/ C u 

C 12 

c 13 

V2C U 

y/2C 15 

V^Ciq 

\ 


C 22 

C 23 

V / 2C*24 

V^C'25 

V 2 C 26 




C 33 

V / 2C < 34 

V 2 C 35 

V 2 C 36 



SYM 


2C 44 

V 2 C 45 

V 2 C 4 Q 






2C*55 

V2C 56 


V 





2C*66 

/ 


where the two indices of the coefficients of matrix Equation (3.27) are written according to the 
notation: C\ \ = Cim, C \2 = Cu 22 , C 13 — C1133, C14 = C1123, C\ 5 = C1131, Ci6 = Cm2, etc. 


3.3.1 Elasticity Equations for an Isotropic Solid 

For an isotropic solid, the stiffness matrix is written: 


Cijkl 3 “ 


where <5^ is the Kronecker delta: 

Sij 


0 if i ± j 
1 if i = j 


(3.28) 


(3.29) 


Lame coefficients are: 


A 


Eu 

(1 + z/)(l — 2u) 


^ 2(1 + v) 

where n is the shear modulus, E is Young’s modulus and v is Poisson’s ratio. 
Substituting , Equations (3.25) and (3.28), into Equation (3.24), we obtain: 

&ij Aefcfcdy + 2/iCjj 


where: 


efcfc = tr(eij) = V • d 


(3.30) 

(3.31) 


(3.32) 

(3.33) 
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Expanding Equation (3.28) into the [6 x 6 ] matrix form yields: 


/ A + 2/z A A 0 0 0 A 

A + 2 fi A 0 0 0 

A + 2 /x 0 0 0 

SYM 2/i 0 0 

2/j 0 

v ; 


(3.34) 


3.3.2 Elasticity Equations for an Ideal Fluid 

For an ideal fluid (nonviscous), the shear stress /i is zero, and the stiffness matrix, therefore, 
reduces to: 


/aaaoooA 

A A A 0 0 0 
A A A 0 0 0 
0 0 0 0 0 0 
0 0 0 0 0 0 
\ 0 0 0000 / 


(3.35) 


Moreover, considering Equation (3.32), we obtain for this specific case: 

(Jij = A tr{eij)5ij (3.36) 

where A is the coefficient of stiffness for an ideal fluid Equation (3.35). Employing Equation 
(3.33) yields 

a = \{V X • d)I 

= -pi (3.37) 

thus verifying 

P = -A(V* • d) 

for an ideal fluid. If we let A be the stiffness parameter we 
in Equation (3.15). 

3.3.3 Dynamic Equations for a Pentamode Material 

In this section we apply the former equations and variables to the specific case of a PM material. 
As mentioned before, a PM has a stiffness matrix which can be written as a [6 x 6] matrix that has 
rank 1 (or only one nonzero eigenvalue). Therefore, considering Equation (3.27), the nonzero 
eigenvalue is the trace of the stiffness matrix C for a PM material becomes: 

eigenvalue = tr( C) = C\ i + C 22 + 6^33 + 2 (C *44 + C 55 + Cq§) (3.39) 


(3.38) 

obtain the fundamental relation stated 
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Moreover, we employ the definition of the stiffness matrix stated in [ 8 ] as the canonical expres¬ 
sion for a PM material: 


C = /CQ ® Q as a A th order tensor 
= /CQQ t as a [6 x 6] matrix 


(3.40) 


where /C has the dimensions of stiffness and is a function of the radial coordinate for the spher¬ 
ical case. For convenience, we rewrite the [3 x 3] tensor matrix Q —» Q, as a single [6 x 1] 
vector defined as: 


Q = 


( Q n \ 

Q 22 
Q 33 
V2Q23 
V2Q31 

V V^Q 12 ) 


(3.41) 


For example, an orthotropic PM, has Q and C given by: 


Q = 


(yen \ 

\[C~22 

y C33 

0 
0 
0 


c = /c 


( c n 

yciiC'22 
yc^c 33 

0 
0 

v 0 


where Q is normalized, so that tr(C) = tr{fC/QC) 1 ) = C\ \ 4 - C 22 + C- 


y C11C22 

C 22 

0 

0 

0 

T\ 


yc 11C33 
y C22C33 

C 33 

0 

0 

0 


0 0 
0 0 
0 0 
0 0 
0 0 
0 0 


33- 


0 \ 
0 
0 
0 
0 
0 / 


(3.42) 


Inserting the canonical expression Equation (3.40) into the stress strain relationship Equation 
(3.24) where the stress are written as [6 x 1] vectors, we obtain: 


a = JCQQ T e 
= {1 CQ t c)Q 


(3.43) 


or 


a = /C[Q 7 e]Q 


(3.44) 


where ICQ T e = (“ scalar ”) in Equation (3.43) and Q and e are [3 x 3] matrices in Equation 
(3.44). We recognize Hooke’s Law applied to a PM material, in the final line of Equation (3.43). 
This expression and the displacement Equation (3.4), incorporating the inertia tensor [ 8 ], here 
rewritten for convenience: 

d 2 d 



constitute the dynamic equations for a PM. 

We note that the last line of Equation (3.43) reduces to Equation (3.37) if Q = I. Moreover, 
in Equation (3.45), g is an inertia tensor, which makes this equation a more general expression, 
applicable to a PM-IC metafluid. 

We define the pseudo-pressure p [8]: 

p = —K, [Q T e] 

= —/C[Qe] (3.46) 

Therefore, the stress tensor Equation (3.43) can be rewritten as: 

a = -pQ (3.47) 

where we note that the tensor Q is symmetric, therefore, Q = Q 7 . Calculating tr(Qe), we 
obtain: 


tr(Qe) — Q\€\ + Q^2 + Q 3 e 3 + 2(^464 + 2Q565 + 2 Q^6 q 
= Q : e 


(3.48) 


Also, observing that, by definition, e = V x d, we obtain the constitutive relationship [8]: 

P = —/C[Q : V^dj (3.49) 


Additionally, the second derivative of Equation (3.49) yields: 


9 2 P _ r 
dt 2 


_ _ fd 2 d\ 

Q: Vx \w) 


(3.50) 


Considering Equations (3.45) and (3.47), we have: 


d 2 d -* 

= • C PQ ) 


(3.51) 


thus, Equation (3.50) becomes: 


d 2 p 


Q V, • (pQ) 


(3.52) 


We note, at this point, that: 


V, • Q = 0 


(3.53) 


must be true in order to be able to compare the “pseudo-pressure” wave equation for a PM with 


25 



the standard wave (Equation (3.1)). Fortunately, this condition can always be satisfied, simply 
by scaling Q. 

Using Equation (3.53), we also note: 


V* • (pQ) = Q V x (p) +p(V, • Q) 
QV x p 


which allows us to write: 


d 2 p 
dt 2 


1C Q : V a 


Q Q V x p 


(3.54) 

(3.55) 


Equation (3.55) is the same as Equation (3.2) and governs the pseudo-pressure (defined by 
Equation (3.46)) within a PM-IC cloak. By a suitable choice of material parameters based 
upon a co-ordinate transformation, it is possible to reduce Equation (3.55) to the regular wave 
equation (Equation (3.1)) within the transformed region. We note that Equation (3.55) becomes 
the pure IC “pseudo-pressure” wave equation (Equation (3.17)), where stiffness is isotropic and 
inertia is anisotropic, if Q = I. Similarly, Equation (3.55) reduces to the pure PM “pseudo¬ 
pressure” wave equation, where the stiffness is anisotropic and the inertia is isotropic, if the 
tensor g is forced to be proportional to I. 


3.3.4 Pure Pentamode Material 

A pure Pentamode material is characterized by an isotropic mass density. Let: 

K=J (3.56) 

and 

g = JQ(V -1 ) 2 Q (3.57) 

which for isotropic density in a pure PM, we require 

6 = JQiV- 1 ) 2 ® = $(*)I (3.58) 

and Equation (3.55) becomes the “pseudo-pressure” wave equation of the material [8]. 

We define h(x) as a scalar function of x : 

K x ) = y® 

such that Va, • (hV) = 0, and assign: 

Q = hV 


(3.59) 

(3.60) 


Consequently, 


g = Jh 2 


(3.61) 
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and the solution for h is obtained if we recall that F = V/ R/ and consider the special case of 
non-rotational deformation (R L = 0). Under this particular assumption, F = F 7 , and hence 
F = V [8]. Moreover, if we recall the identity of Equation (3.13) and recall the constraint 
Equation (3.53) we see that for F = V in Equation (3.60) the constraint is met if we let: 


Q = J- 1 F (3.62) 

Equation (3.58) becomes: 

6 = JQ(V 1 ) 2 Q 
= J(J' 1 F)F" 2 J' 1 F 

= J- L I (3.63) 


Consequently, the stiffness constant and the isotropic inertia elements become, respectively: 

/C = JJC o (3.64) 

q = poJ- 1 ! (3.65) 

where /Co and p 0 are the outer fluid constants bulk modulus and density, as previously for the 
IC metafluids case. 

Considering Equation (3.62) and substituting Equations (3.64) and (3.65) in Equation (3.55), 
we obtain the wave equation for a pure PM material without coordinate rotation: 

— F : V*[F V x p] (3.66) 

or 

= F r : V x [F T V x p] (3.67) 

3.3.5 Boundary Conditions for a PM-IC metafluid 


The boundary conditions for a PM material follow the same lines as for the IC metafluid, re¬ 
garding the continuity of normal stress and displacement. 

The pseudo-pressure p must match the pressure field of the exterior fluid, once the transfor¬ 
mation is performed. Therefore, we impose that the normal component of the stress tensor, 
Equation (3.47), must satisfy the condition, at r = b [8]: 

[nQp] = 0 (3.68) 

Regarding the continuity of the displacement we consider the normal component of Equa- 
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tion(3.51), noting V* • Q, and imposing, at r = b: 

[h • g-'QVvp] = 0 (3.69) 

At the inner surface of the cloak, the conditions of continuity in the pressure field and displace¬ 
ment are also imposed. 

3.4 SPHERICAL GEOMETRY 

In this section, we address the application of the acoustic cloaking theory to the specific case of 
a spherical geometry. 

We fix polar and azimuthal angles. Hence, only the radial coordinate r need be considered, and 
orthogonal coordinates 6 and 0 are unchanged. 

We consider the inverse transformation [8]: 

X = f(r)x (3.70) 


where x = x/r, r = |x|. 

^ ^ 

The function /(r) defines a transformation of coordinates from x to X, and we force it to be 
monotonic in the interval a < r < b and be such that: 


f(r) 


6 if r < a 

b if r — b 

r if r > b 


(3.71) 


We define radial and tangential matrices in Cartesian coordinates (I r and I ), such that: 


I = I r + I± 


where 


I r = x <8) x 

= XX 


and 


I_i_ = I — x ® x 
= I — xx T 


(3.72) 


(3.73) 


(3.74) 


Moreover, for the general case of a spherically symmetric acoustic metafluid, the symmetric 
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stress tensor Q is of the form [ 8 ]: 


Q = cu(r)(I r + 7 (r )I_l) (3.75) 

where u> and 7 are scaling functions depending on r. 

Cloak parameters IC, Q and o are specifically defined for IC and PM in the following paragraphs. 


3.4.1 Inertial Cloak Applied to the Spherical Geometry 

In order to define the IC parameters applied to the spherical geometry instance, we determine 
parameters given by Equations (3.18) and (3.19). 


Considering Equations (3.5) and (3.6), and Equations (3.73) and (3.74), the deformation tensor 
F becomes: 

F iJ — 777 \Ir + r, (3.76) 

f (r) f(r) 

where f'(r) = Furthermore, the determinant J of F is 


J 


f(r)f 2 (r) 


(3.77) 


We have shown that the stress tensor of an IC metafluid is isotropic, hence Q 
Va; • Q = 0 is immediately satisfied. 

Using Equation (3.77), Equations (3.18) and (3.19) can be written as: 

„2 


/c = 


r 


f ( r ) f 2 ( r 


-IC 


0 


r 

q = P0777 \ f2( \ (y 
f ( r)f 2 (r ) 


—1 \2 


Using FF T = Vl and Equation (3.76): 


I, and condition 


(3.78) 

(3.79) 


(V - 1 ) 2 


[f\r)]\ + 



L r 


(3.80) 


The anisotropic inertia tensor is of the form g = p r I r + p±l± [ 8 ]. Applying Equation (3.80) 
into Equation (3.81), we obtain: 


Q — P 0 



/» C + 



(3.81) 
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We note that it is possible to define /C and g as functions of the radial and tangential components 
of the speed of sound within the cloak. 

The wave velocity in the exterior fluid is defined by Cq = /Co/p 0 (m/s), as mentioned before. 
Dividing Equations (3.78) and (3.81), we find the wave velocity within the IC metafluid: 


ClC — 1 — ^Q)Po 


-1 


I r + 


U/'(' r )J r L f(r) 


(3.82) 


where the radial and tangential components of the wave velocity: c r = \J /C/ p r and c_l = 
VWpI can be identified from Equation (3.82): 

c, = (3.83) 

/ M 

and 

C -L = ~TT~\ (3-84) 

f(r) 

If c r = c± the cloak density becomes isotropic [8]. 

At the outer surface of the cloak r = b, the unit vectors obey: h = N — x — e r and the 
pressure field must be continuous, [p] = 0. This condition is satisfied if we consider Equation 
(3.23) at r = b. Additionally, the IC parameters must match the exterior fluid constant bulk 
modulus and mass density: JC = JCq and p = p 0 , which means that conservation of mass must 
be obeyed. We consider the radial component of Equations (3.76) and (3.81) at r = b: 


F ,.(&) = 


1 

m 


(3.85) 


Qrip) = Po 



(3.86) 


where we have used Equation (3.71): f(b) = b. Therefore, we are able to confirm, using the 
equality given by Equation (3.20) ff(p 0 p _1 ) = iVF at r = b, which ensures that the mass 
density is preserved. 


3.4.2 Pentamode Applied to the Spherical Geometry 

Because of spherical symmetry and the lack of rotation in the coordinate transformation, the 
PM-IC metafluids have transverse isotropic (TI) symmetry [20]. A regular solid with TI sym¬ 
metry comprises five independent elastic moduli where X:> is an axis of symmetry, hence the 
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solid is isotropic in x\ and x -2 directions [ 8 ]. The stiffness matrix is defined as follows: 


( C n 

C12 

C13 

0 

0 

0 

\ 

Cl 2 

Cn 

C13 

0 

0 

0 


Cy, 

c 13 

C33 

0 

0 

0 


0 

0 

0 

C44 

0 

0 


0 

0 

0 

0 

C44 

0 


V 0 

0 

0 

0 

0 

|(Cn-C 12 ) 

/ 


(3.87) 


A PM-IC metafluid has only two independent elastic moduli, has no shear stress, and has only 
one proper solution. In order to obey these conditions, we let C 44 = 0 and Cn = Ci 2 . Cijki 
defined in Equation (3.87) becomes: 


/ 

C n 

C n 

C13 

0 

0 

0 

\ 


C n 

Cn 

C13 

0 

0 

0 



c 13 

c 13 

C33 

0 

0 

0 



0 

0 

0 

0 

0 

0 



0 

0 

0 

0 

0 

0 


V 

0 

0 

0 

0 

0 

0 

/ 


(3.88) 


and C 4 3 — \J C'i iC' 33 . 


The definition of C t jki in Equation (3.88) satisfies the definition of Cijki given by Equation 
(3.42). We let Cf 3 = K r , C\ x = K± [ 8 ]. 


Considering Equation (3.75), for the specific instance of a pure PM material u and 7 are defined 
as follows [ 8 ]: 



(3.89) 


l r 


and 

/ \ f'(r) 

7[r) = —^~ r 

and the condition V* • Q = 0 will be satisfied. 
Therefore,we may write Equation (3.40) as follows: 


(3.90) 


C = /C 
= K 


ui(r)I r + tu(r)7(r)Ij_ tu(r)I r + cu(r)7(r)Ij_]' 
u 2 (r) I r + cu(r) 2 7 2 (r)I ± j 


(3.91) 


where the symmetric property between tensors I, and I have been applied. Radial and tangen- 


31 



tial anisotropic stiffness components, for a pure PM, are given by [8]: 


K r 


/Cce 2 (r) 


1 

f'(r) 


7 ( 01 2 

- r . 


and 


K_ l = K,u 2 (r) 7 2 (r) 

= /(0 


The isotropic inertia tensor is given by [8]: 


Q = Po/V) 



We note that the relationship between the IC and the PM parameters: 


(3.92) 


(3.93) 


(3.94) 


{ IC, p r ■ P± } —> {1/p, 1//C r , 1//C_L }. 


The wave velocity expression for a PM material is also found, as for IC, by performing the 
division between the stiffness components and the isotropic density vector. The radial and 
tangential components of the wave velocity for a pure PM material are: 


c r = 


1 

f(r) 


(3.95) 


and 

c ± = T7“T (3-96) 

f(r) 

Equations (3.95) and (3.96) match Equations (3.83) and (3.84), which confirms that the wave 
velocity parameters do not depend upon the metafluid we are considering, as anticipated in the 
acoustic cloaking theory section [8]. 


3.4.3 CLOAKING THE SPHERICAL SHELL 

In this section we reduce the “pseudo-pressure” wave equation to the regular wave equation, in 
order to obtain the solution for the acoustic pressure field within the cloak (pseudo-pressure). 
The net effect of the transformation of coordinates results in a transformed space R = f(r), 
where /(r) is defined in Equation (3.71). The radius of the shell has been shrunk from r = a to 
R — S, where 0 < <5 < a. 
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We consider the pseudo-pressure wave Equation (3.55), here rewritten for convenience: 


d 2 p 
d t 2 


1C Q : V a 


' Q 


(3.97) 


We want to define the right hand side of Equation (3.55) in spherical coordinates, and in the 
transformed coordinate system R, where the solution for the pseudo-pressure field can be writ¬ 
ten in the regular form of the wave equation solution. 


Considering the unit vectors in spherical coordinates: 


s r 


e e 


&(j> 


sin 0cos 0 
sin 0sin 0 
cos 9 

cos 9 cos 0 
cos 9sin 0 

— sin 6 

— sin 0 
COS 0 

0 


(3.98) 


The Laplacian operator applied to p yields: 


A dp „ 1 dp „ l dp 
Va jP = 6 r — + e g -— + —T —- — 

or r 09 r sin 9 dcp 


(3.99) 


The tensor Q was defined in Equation (3.75) and oj{r) and 7 (r) are given by Equations (3.89) 
and (3.90). 


We multiply Q S7 x p and obtain: 

dv 1 dv 1 dv 

QVa -P = u(r)-^-e r +u(r)j(r)-^e e + u(r)j(r) —— (3.100) 

or r d9 r sin 9 dtp 

where we note e r = x, I r = e r ej and I = I — e r ej. Therefore, I r e r = e r , Ig£ r = 0 and 
Ipe r = 0, because the unit vectors e r , eg, e$ are mutually orthogonal. 


Using Equation (3.94), we find: 


q 1 QS7 x p 


1 dp „ 1 dp „ 1 dp „ 

Pof'(r ) Or 6r + p 0 f (r) <90 60 + p 0 f(r) sin0 <90^ 


(3.101) 


The product /CQ is calculated using Equations (3.77) and again, Equation (3.75). We finally 
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obtain: 


/CQ : V. o x QV x p 


/C 0 r 1 <9 / 1 9p\ 2 /9p\ 

p 0 l/X'r) 9r ^/ ( r ) dr) /(r)/'(r) \9r/ 

1 <9 / .gg\ 1 /9p\l 

/ 2 (r) sin 9 99 v 99/ / 2 (r) sin 2 9 \90// 


(3.102) 


We employ the transformational function /? = /(r) into the right hand side of Equation (3.102) 
and substitute the result into Equation (3.97): 

8 2 p JC 0 c 8 dp \ 2/(9 p\ 1 8 / 9p\ 1 / 9p \ t 

8t 2 p 0 ^ <99? <97?/ + 7? V<97?/ + 7? 2 sin 9 <96* V Sm <96*/ + 7? 2 sin 2 9 V90/ / 


where we have: 


or 


0 



dr 8 18 

(97? dr f (r) dr 


8 2 p _JC 0 f 1 9 / 2 dp\ 1 9 / . (9p\ 1 /9p\) 

8t 2 p 0 17? 2 97? V 97?/ 1 7? 2 sin9 99 V sm 997 7? 2 sin 2 9 V907 J 


(3.103) 


(3.104) 


We recognize the Laplacian operator at the right hand side of Equation (3.104), and write this 
equation as: 


8 2 p 

W 


ICo 

Po 


(3.105) 


Equation (3.105) is the standard wave equation for the pseudo-pressure, within the transformed 
region. The solution of this equation assumes the same form as the solution of the pressure 
field at the exterior cloak, because it is the standard solution for the wave equation in spherical 
coordinates. For a thin spherical shell, this result was determined in Chapter 2, in Equation 
(2.34). For this specific case, yields: 


p(r) = ^2 [®njn(kff(r)) + p n h n (k f f(r))\P n (cos0) (3.106) 

?i=0 

where a n and 3 n denote the pressure amplitude scaling modal coefficients. 

The function f(r) must be defined in order to satisfy the boundary conditions already depicted 
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in Equation (3.71). Additionally, we will impose the condition 


/>) = 1 


(3.107) 


This condition is not a fundamental requirement in the definition of an acoustics cloak. The 
fundamental requirements lie in the condition that f(b) = b and that f(r) be monotonic, while 
continuity of the normal velocity and stress at r = b is imposed, forcing impedance matching 
of the external fluid and acoustic cloak at r = b, however, employing Equation (3.107), we are 
introducing an additional requirement at the outer boundary of the cloak, ensuring a smooth 
transition in density between the two media. 

In order to impose the boundary conditions stated in Equations (3.71) and (3.107), the function 
f(r) must be quadratic. 


Let 


f(r) = 5 + A(r — a) + B{r — a) 2 


applying Equation (3.107) we determine 


(3.108) 


1 - A 
2(6 — a) 


Using this result and applying condition f(b) = b we determine 


A 



- 1 


Substituting these results into Equation (3.108), we obtain: 


f(r) = 5 + (6 + a 


25) 


~r — a~ 

.b — a. 


+ (S 


a) 


~r — ai 2 

.b — a. 


(3.109) 


The coefficients a n and f3 n in Equation (3.106) must also be determined. The boundary condi¬ 
tions of the system are, again, employed to perform this evaluation. As mentioned before, the 
variation of the total acoustic pressure field and the normal particle velocity have to match at 

the outer and inner surfaces of the acoustic cloak. We obtain the normal component coefficient 

—^ 

W n of the displacement vector d considering the radial component of Equation (3.51): 


-co 2 W(r) 


1 dp 
pof(r) dr 


and evaluate W at r = a + and at r — b : 


W (a + ) 


1 dp 
poU 2 f'(a) dr 


(3.110) 


(3.111) 
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and 


W{b~) 


1 dp 
p 0 u 2 f'(b) dr 


(3.112) 


Moreover, considering the boundary condition/(6) = b, and taking the amplitude coefficients 
of Equations (2.19) and (3.106), we may perform the equalities: 

Pn S = ot n j n (kr) + /3 n h n (kr) 

or 

(2 n + 1 )i n j n (kr) + R n h^\kr) = a n j n (kr ) + )3 n h n (kr ) (3.113) 

Therefore, a n = (2 n + 1 )i n and /3 n = R n - Equation (3.106) becomes: 

OO 

p(r) = + 1 )* n in(^//(»’)) + Rnhn(kff(r))]P n (cos9) (3.114) 

?i=0 

and the radial displacement: 

1 X 

W(r) = -J^[(2n + l)i n j' n {k f f(r )) + R n h' n (k f f(r))\P n (cos9) (3.115) 

P o w J \ r > n—n 


In order to ensure the continuity of the normal displacement at r = a, we consider W (a) and 
evaluate the scaling amplitude coefficient R n at this position: 

Rnh' n (k f f(a )) = uz a W n - (2 n + l)i n j' n (k f f(a)) (3.116) 

where we defined the specific acoustical impedance z a coming from the relations: = 

i _ i 

UJpQCf UlZa 


Also, considering the equation of the total pressure field for the thin shell Equation (2.19) and 
the definition of the specific acoustical impedance, Equation (2.32), we can state: 

-iuz n W n = (2 n + 1 )i n jn(k f f(a)) + R n h n (k f f(a )) (3.117) 

Thus, 

R n h n (k f f(a )) = -iuz n W n - (2 n + l)i n j n (k f f(a )) (3.118) 


We multiply iz n by Equation (3.116) and Z a by Equation (3.118) and sum these two results, in 
order to eliminate W n , and obtain the definition of R n : 


(2 , lV -n iZnj' n (kff(a )) + z a j n (k f f(a )) 

iz n h' n (k f f(a )) + z a h n (k f f(a )) 


(3.119) 


36 



where we note that as 4 —» 00, Equation (3.119) becomes the expression Equation (2.11), 
which is the scaling coefficient for an acoustically rigid boundary condition. Similarly, as z n —>■ 
0, Equation (3.119) becomes Equation (2.16), which is the expression for the scaling coefficient 
of the “pressure release” boundary condition. 

The pseudo-pressure is, therefore, defined by substituting Equation (3.119) into Equation (3.114) 


p(r) = ^[(2n + 1 )i n 

n =0 


3n(k f f(r )) 


iZnj'n(k f f(a)) + Zaj n (kff(a)) 

iz n h' n (k f f(a)) + z a h n (k f f(a)) n[ fJ[ >J 


P n ( cos6>) 

(3.120) 
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CHAPTER 4: 
DATA ANALYSIS 


In this chapter we develop an ensemble of experiments performed in MATLAB in order to 
evaluate the effectiveness of the acoustic cloak applied to the thin spherical shell. 

4.1 SCATTERING ENERGY 

The most interesting results, for our purposes, are related to the evaluation of the amplitude of 
the scattered pressure field in the far field. We employ form functions [10] and [21], which can 
represent the scattered pressure amplitude as a function of a single direction (a particular value 
of the polar angle 6), or a scattering profile in all directions. 

The form functions are defined by: 


/c 


lim 

r—>oo 


2r\p s (r,9)\ 

aPi 


(4.1) 


where | ps(r, 6 ) | denotes the scattered amplitude and P, is the pressure amplitude of the incident 
pressure, as defined in Chapter 2. 


4.1.1 Scattering Energy Coefficients 

The scattering energy coefficients [11] are scalar quantities that provide an evaluation similar 
to a figure of merit, that register the total scattered energy over the entire surface of the sphere. 
This result is very interesting in the sense that it allows us to evaluate and compare the scattered 
energy before and after the application of the cloaking material. 

The acoustic scattering coefficients are defined by: 

4 °° 

< 2 » + ') 

v ' n=0 

where ka is the product of the wave number and the radius of the spherical object, and b mn 
and a mn are the scaling coefficients of the scattered pressure field amplitude, determined in 
the previous chapters for each particular considered problem: hard and soft spheres and the 
spherical shell without and with a cloak. 

The definition of the scattering coefficients can be extended from the spherical geometry appli¬ 
cation to other objects of finite shape. A closer analysis of this parameter shows that it provides 
a measure of the energy scattered by the projection of the surface area on to a plane perpendic¬ 
ular the direction of propagation of the incident plane wave. Consequently, a c is a very good 
estimate of the forward scattered energy. 


(4.2) 
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4.1.2 Model System Parameters 

We evaluate our model spherical obstacle/medium system by considering the following condi¬ 
tions: the surrounding medium has the acoustic wave speed and density of the sea water and the 
spherical obstacle is made of steel, of radius lm. For the shell case, we assume a thin spherical 
shell, with thickness 0.05m and mid-surface radius a = lm. Moreover, we consider the ideal 
case of negligible sound pressure waves within the structure of the shell, and assume that there 
is no variation of the spherical shape due to torsional motions or variations of the mid-surface 
radius. The outer surface of the cloak has radius 6 = 2. These assumptions consist of an ideal 
formulation; nonetheless, they are convenient in the sense that we want to evaluate the perfor¬ 
mance of the acoustic cloak. The amplitude of the incident plane pressure wave is IPa. We 
will consider 2 types of cloaks: pure Inertial Cloak (IC) and a pure Pentamode (PM). 

We perform the analysis of the system and determine the scattered pressure patterns and other 
parameters of interest considering the low frequencies of 150 Hz, 1500 Hz and 5000 Hz. Ta¬ 
ble 4.1 shows the natural frequencies of the submerged thin spherical shell, with no cloak. 
These values were equating using the characteristic equation for the loaded spherical shell [10]. 

Table 4.1: Natural frequencies of submerged thin spherical shell 


n th Normal Mode 

ka 

f w (Hz) 

ka 

f (2) (Hz) 

0 

— 

— 

5.36 

1280.17 

1 

0 

0 

6.72 

1603.74 

2 

1.82 

434.70 

9.42 

2249.4 

3 

2.19 

523.81 

12.65 

3020.1 

4 

2.45 

584.82 

16.03 

3826.1 

5 

2.62 

625.80 

19.46 

4644.7 

6 

2.74 

653.30 

22.91 

5469.1 

7 

2.82 

672.92 

26.38 

6296.5 

8 

2.89 

687.73 

29.85 

7125.9 

9 

2.93 

699.36 

33.33 

7956.2 






oo 

3.36 

800 

oo 

oo 


The natural frequencies are related to 12 by 

c 12 
2na ( c/cp ) 


8.3511 


The breathing mode (n = 0) has only one resonant frequency, and the lower resonant frequency 
for n — 1 is zero (first root of the characteristic equation equals zero). Higher modes have two 
resonant frequencies. Modes higher than n — 5 present ka > 20.9, that is, the frequency is 
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higher than 5000 Hz and are not reflected in our scattering pattern analysis. Nonetheless, the 
lower-frequency modes are the ones with more radial motion, hence, that better couple with the 
acoustic field. 


4.1.3 Acoustics Scattering Pattern 

We begin our analysis of data by observing the scattered energy when there is no application 
of a cloak. Figure 4.1, Figure 4.2 and Figure 4.3 depict the scattering pressure pattern from a 
homogeneous sphere in the two limiting impedance conditions, rigid and pressure release, as 
well as the scattering pattern for the thin shell. The incident pressure wave is coming from the 
180° direction. Therefore, 180° represents the direction of the backscattered energy, and the 0° 
represents the forward energy direction. 


Scattered Pressure Field. Far Field. Thickness of the Shell 0.05 
11 Mode(s), summed from n = 0 to n =10 
Frequency 150 Hz 




Figure 4.1: Pressure scattered amplitude pattern from a spherical shell and a homogeneous 
sphere - rigid and pressure release boundary conditions. Frequency 150 Hz, ka = 0.268. 
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Scattered Pressure Field. Far Field. Thickness of the Shell 0.05 
11 Mode(s), summed from n = 0 to n =10 
Frequency 1500 Hz 

90 

20 



- - - Rigid 

Press.Release 
Sph. Shell 


0 


Figure 4.2: Pressure scattered amplitude pattern from a spherical shell and a homogeneous 
sphere - rigid and pressure release boundary conditions. Frequency 1500 Hz, ka = 2.68. 


Scattered Pressure Field. Far Field. Thickness of the Shell 0.05 
31 Mode(s), summed from n = 0 to n =30 
Frequency 5000Hz 

90 


25 



- Rigid 

Press.Release 
Sph. Shell 


0 


Figure 4.3: Pressure scattered amplitude pattern from a spherical shell and a homogeneous 
sphere - rigid and pressure release boundary conditions. Frequency 5000 Hz, ka = 20.9. 
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We observe that, at this low frequency, the “pressure release” scattered patterns occur in all 
directions, the rigid sphere presents almost only backscattered energy and the spherical shell 
scatters in the azimuthal directions, parallel to the incident plane (assuming a pinched balloon 
shape), suggesting back and forward vibrations. 

As the frequency increases, so does the directivity, narrowing the main lobe and increasing 
the amplitude of the side lobes, as expected, for all the three cases: acoustically hard and soft 
spheres, as well as the spherical shell. We also see that the pressure scattering amplitude of the 
shell increases significantly at 1500 Hz, which lets us infer that we are close to the fundamental 
frequency, or one of its harmonics, for this specific case. 

4.2 APPLICATION OF THE CLOAK 

4.2.1 Effectiveness of the Acoustic Cloak 

In this section we apply the cloak to our spherical shell system. The application of the transfor¬ 
mational acoustics method has the net effect of making the spherical obstacle appear smaller to 
the incident plane wave, as explained in Chapter 3. Figure 4.4, Figure 4.5 and Figure 4.6 depict 
the same scattering patterns of the shell as seen in the previous section, but with a cloak applied, 
for transformed radii of <5 = 0.1a, S = 0.5a, and S = 0.95a. The effectiveness of the acoustic 
cloak for each case is given. 


Scattered Pressure Field. Thickness of the Shell 0.05 
11 Mode(s), summed from n = 0 to n =10 
Frequency 150 Hz 

90 


5 



5 = 0.5 
5 = 0.95 
no cloak 


0 


Figure 4.4: Pressure scattered amplitude pattern from a spherical shell without and with cloak. 
Net effect of the application of the cloak as the transformed radius is reduce to <5 = 0.1a, 
5 = 0.5a or <5 = 0.95a. Frequency 150 Hz, ka = 0.268. 
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Scattered Pressure Field. Thickness of the Shell 0.05 
11 Mode(s), summed from n = 0 to n =10 
Frequency 1500 Hz 

90 

20 



- 5 = 0.1 

- 5=0.5 

- 5=0.95 

- no cloak 


0 


Figure 4.5: Pressure scattered amplitude pattern from a spherical shell without and with cloak. 
Net effect of the application of the cloak as the transformed radius is reduce to S = 0.1a, 
S = 0.5a or 5 = 0.95a. Frequency 1500 Hz, ka = 2.68. 


Scattered Pressure Field. Thickness of the Shell 0.05 
31 Mode(s), summed from n = 0 to n =30 
Frequency 5000 Hz 


90 



Figure 4.6: Pressure scattered amplitude pattern from a spherical shell without and with cloak. 
Net effect of the application of the cloak as the transformed radius is reduced to <5 = 0.1a, 
S = 0.5a or 5 = 0.95a. Frequency 5000Hz, ka = 20.9. 
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We can observe that the pressure scattering amplitudes of the shell for the case of the frequencies 
150Hz and 1500Hz have been considerably reduced with the application of the cloak. Even for 
a reduction of the transformed radius to only 0.95a, the application of the acoustic cloak causes 
the pressure scattered amplitudes to be reduced by a ratio of roughly 1/40 and 1/3, respectively, 
which represents a very good effectiveness of the acoustic cloak at these lower frequencies. 
Nonetheless, we also observe that the decrease in the pressure scattered amplitude does not 
happen for a frequency of 5000 Hz. At this frequency, the application of the cloak seems to be 
effective only if the radius of the shell is effectively reduced to 0.1a. Unfortunately, decreasing 
the effective radius causes the physical parameters (mass density and stiffness components) to 
reach values that are outside the bounds of physical implementation. It is important to note, 
however, that the application of an acoustic cloak that simulates a reduction of the spherical 
shell radius to an half of its original value (6 = 0.5a) causes the forward scattered pressure 
to be higher than the original (with no cloak). On the other hand, the backscattered pressure is 
significantly reduced, which can provide interesting operational applications. Table 4.2 presents 
the scattering coefficients as defined by Equation (4.2) for the homogeneous sphere at the two 
impedance limiting boundary conditions and for the spherical shell, at various transformed radii 
150Hz, 1500 Hz, 5000 Hz frequencies. 

Table 4.2: Scattering coefficients for the homogeneous sphere - rigid and “pressure release” boundary conditions, 
and for the spherical shell - without and with cloak at different radius lengths 


Scattering Coefficients (a) 

Frequency (Hz) 

150 

1500 

5000 

Rigid Sphere a 

0.0851 

1.4035 

1.7460 

“Press-Release” Sphere a 

3.6384 

2.5264 

2.2501 

Shell a 

6.1842 

4.6486 

2.0010 

Shell 00.95 

0.0224 

1.7401 

2.1040 

Shell (To .05 

2.67 x 10“ 4 

0.1362 

0.5180 

Shell o 'o o! 

1.75 x 10“ 7 

0.0128 

0.0106 


We observe that the scattering coefficients of the shell, with no cloak, are in general higher than 
the homogeneous sphere. Furthermore, the scattering coefficients’ values drop significantly as 
the radius of the shell is decreased, at 150 Hz and 1500 Hz. At 5000 Hz we verify the need of 
further reducing the radius of the shell. These results are consistent with figures 4.1 to 4.3 and 
4.4 to 4.6. 

4.2.2 Normal Modes 

We analyze the scattered pressure amplitude of the shell with and without the cloak in terms 
of the fluid-loaded normal modes in order to evaluate which modes radiate more energy and 
the effectiveness of the cloak for each mode. We note that for the specific case of a spherical 
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geometry, the normal modes are independent of each other (which makes possible this exami¬ 
nation). Figure 4.7, Figure 4.8 and Figure 4.9 depict the scattered pressure amplitudes in terms 
of the normal modes for the spherical shell case, without and with cloak at the transformed 
radius lengths 5 = 0.95, <5 = 0.05 and 5 = 0.01. We note that these modal scattering values are 
complex. Therefore, when the scattered modal fields are summed cancellation occurs. 


Scattered Pressure Field. Far Field. 
11 Mode(s), from n = 0 to n =10 
Frequency 150 Hz 



Figure 4.7: Pressure scattered amplitude from a spherical shell in terms of normal modes. Fre¬ 
quency 15077,2, ka = 0.268. 
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Scattered Pressure Field. Far Field. 
11 Mode(s), from n = 0 to n =10 
Frequency 1500 Hz 



Figure 4.8: Pressure scattered amplitude from a spherical shell in terms of normal modes. Fre 
quency 1500-ffz, ka = 2.68. 


Scattered Pressure Field. Far Field. 
31 Mode(s), from n = 0 to n =30 
Frequency 5000 Hz 



Figure 4.9: Pressure scattered amplitude from a spherical shell in terms of normal modes. Fre 
quency 5000Hz, ka = 20.9. 
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In Figure 4.7 and Figure 4.8 we can see that the scattered energy is concentrated roughly in 
the first 2 and 6 lower modes, at 15077c; and 150077c; frequencies, respectively. This is the 
qualitatively expected result for these lower frequencies. 

At 500077c;, for <5 = 0.1a and <5 = 0.5a, the scattered energy is spread out to roughly the first 
10 modes, for the uncloaked and cloaked shells. These results are consistent with what was 
observed for 15077c; and 15007 Iz. However, for <5 = 0.95a, a number of unexpected results 
are observed. For one, we see that the first 24 modes are excited, which is more than twice 
the number that was excited for the uncloaked case, and for the cloaked cases where 5 = 0.1a 
and d = 0.5a. Consequently, one might expect greater scattering than that for the other cases, 
based upon the modal amplitudes, yet this is not observed in Figure 4.5. As already noted, the 
modal amplitudes are complex quantities and add in and out of phase. An understanding of 
these results is not obvious and requires more research. 

4.2.3 Backscattered Energy 

We investigate the behavior of the cloak in the backscattered direction (180° direction), consid¬ 
ering a range of frequencies from 0Hz to 500077 z. Figure 4.10 illustrates the variation of the 
backscattered energy in terms of ka, where ka is the dimensionless frequency. 


Comparisson between backscattered amplitudes 



ka 


Figure 4.10: Pressure back scattered amplitude of a homogeneous sphere - rigid and pressure 
release boundary conditions and for a spherical shell, with and without cloak - at 5 = 0.95, 
S = 0.05 and 5 = 0.01 as a function of ka. 


Figure 4.10 allows us to determine the fundamental frequencies and the respective harmonics 
for the conditions we are investigating. It confirms that the application of the cloak causes a 
shift of the fundamental frequencies to higher values. We observe that this shift is inversely 
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proportional to the transformed radius length. In other words, we verify that, in the backscat- 
tered direction, as the transformed radius is decreased, the resonance frequencies are shifted 
to higher frequencies. This effect is possibly due to the fact that as the radius of the shell is 
decreased, the modal frequencies are being evaluated at k5, rather than ka. Therefore, the fre¬ 
quency to line up the spikes must be increased by a factor of (a/5). Moreover, we observe a 
decrease in amplitude proportional to 5. By the same argument, and considering the scaling 
factor 2/ka in Equation (4.1) the amplitude of the cloaked curves must be multiplied by S, in 
order to obtain the non cloaked pressure amplitudes. 

4.2.4 Acoustic Cloak Parameters 

This section considers the cloaking parameters, mass density and stiffness coefficients, within 
the acoustic cloaks. In Figure 4.11 we consider a reduction of the transformed radius of the shell 
to half of its original value. This figure displays the relative magnitude of these parameters 
between the inner (r — a — 1 m) and the outer (r = b = 2m) boundaries of the acoustic 
metafluid. 


Normalized Material Parameters. 
Thickness of the Shell 0.05 
Transformed radius (5) :0.5a 



Figure 4.11: Acoustic cloak parameters. Frequency 150 Hz, ka = 0.268. 


The acoustic parameters are normalized by sea water density and stiffness component’ values. 
We observe that at r = b all the parameters approach unity, providing a smooth transition at the 
outer boundary of the cloak. At r = a we highlight the result mentioned earlier in this work. 
We verify that the radial anisotropic mass density of a pure IC metafluid reaches a maximum 
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as the inner surface of the cloak is approached. Similarly, the radial anisotropic stiffness com¬ 
ponents of a pure PM material tend to zero at the inner surface of the cloak. Table 4.3 and 
Table 4.4 depict the numerical values for the density and stiffness modulus, regarding different 
transformed radius lengths. The values were determined using Equations (3.78), (3.81), (3.92), 
(3.93), (3.94) and (3.109). 


Table 4.3: Inertial cloak parameters 


Transformed Radius 

INERTIAL CLOAK 

Prip^) 

Pr(b) 

P±(a) 

P±(P) 


Kiso{b) 

0.95 

1.2 

1 

0.9 

1 

l 

1 

0.5 

8 

1 

0.5 

1 

2 

1 

0.1 

280 

1 

0.4 

1 

35.7 

1 


Table 4.4: Pentamode parameters 


Transformed Radius 

PENTAMODE 

lC r (a) 

K r {b) 

JC±(a) 

M&) 

Pisoifl) 

Pisoip ) 

0.95 

0.8 

1 

1.1 

i 

1 

l 

0.5 

0.1 

1 

2 

i 

0.5 

l 

0.1 

0.004 

1 

2.8 

i 

0.03 

l 


We highlight the value obtained for the radial density of an IC at the inner boundary of the 
acoustic cloak r — a (Table 4.3), p r (a) = 280 and the value obtained for the radial stiffness 
coefficient, also at r = a, JC r (a) = 0.004 (Table 4.4). These values reinforce the asymptotic 
behavior of the acoustic cloak parameters showed in Figure 4.11 at r = a. There is no physical 
feasibility of a material with a density 280 greater than the density of water, or a material with 
zero stiffness. 
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CHAPTER 5: 
CONCLUSIONS 


In this thesis we were able to analytically solve the scattering problem and use the transforma¬ 
tional acoustics method applied to a spherical geometry to investigate acoustic cloaking. While 
the analytical development is only possible for well defined geometric shapes of the obstacle (or 
the target), such as a sphere, it enabled us to investigate the fundamental theoretical details that 
constrain the final solution as well as to acquire a more precise notion of the variables involved 
and their effect on cloaking performance. 

We began our work in Chapter 2 by determining the incident and scattered pressure amplitudes 
for a homogeneous spherical object with two limiting impedance boundary conditions: rigid 
and “pressure release.” This approach allowed us to state and bound our problem. The scattered 
pressure amplitude from a thin spherical shell was then determined. The scattered pressure 
amplitudes for the different case studies considered are all defined in spherical coordinates, in 
terms of spherical Bessel functions and spherical harmonics expansions. They are scaled by 
the incident pressure amplitude, which allows us to compare the different results and scattered 
patterns obtained. We observed that, at lower frequencies, the spherical shell registers greater 
scattered amplitudes than the homogeneous sphere, as well as a different radiation pattern. As 
the exciting (incident) frequency increases the main lobe beam-width narrows down. Therefore, 
the scattered directivity increases, as well as the side lobes’ amplitude and width. 

Chapter 3 is dedicated to the transformational acoustics development, in order to apply an acous¬ 
tic cloak to the spherical object. The transformation of coordinates between an original coor¬ 
dinate system and a deformed coordinate system is explored and the general acoustic cloak 
theory is developed. The original coordinate system contains a thin spherical shell immersed in 
a metafluid that is anisotropic in density, anisotropic in stiffness or anisotropic in density and 
stiffness. The transformational acoustics method simulates shrinking the radial length of the 
spherical shell. This results in a smaller scattered surface area, thereby, lowering the scattered 
energy. The outer surface of the acoustic cloak remains the same (r = b ) and the normal dis¬ 
placement and particle velocities of the metafluid and the outer fluid (assumed sea water) match 
along this surface. 

Two different theoretical types of acoustic metafluids were examined that varied in their physi¬ 
cal properties of density and stiffness. These metafluids were considered to be ideal, inasmuch 
as the anisotropic properties varied continuously. The IC is easier to explore analytically and to 
model in a computer simulation . Using the definition of the pressure field in terms of the stiff¬ 
ness and divergence of the displacement, and the linearized Euler’s equation, we were able to 
define the IC cloak parameters, anisotropic density and isotropic stiffness, as well as determine 
the limitations of the physical implementation of this metafluid when the boundary conditions 
are imposed. The pentamode materials are designed by considering first the stress strain re¬ 
lationship and determining the stiffness matrix needed to obtain and characterize the pressure 


51 




field (“pseudo-pressure”) within this metafluid. The assumption of a divergence free, non¬ 
unique, tensor Q makes it possible to define the necessary stiffness matrix and consequently 
the PM parameters for an isotropic density with anisotropic stiffnesses. Again, imposition of 
the boundary conditions allowed us to determine the limitations on the physical implementation 
of this type of acoustic cloak. Additionally, assignment of the acoustic cloak parameters of 
density and stiffness for IC and PM showed that the sound speed c (m/s) is the same in either 
case, which allows us to conclude that it would have been possible to define the acoustic cloak 
parameters in terms of the radial and perpendicular components of c. 

In Chapter 4 we developed some computational examples in order to evaluate the effectiveness 
of the acoustic cloak, in terms of the frequency and the transformed radius 5. We explore the 
polar diagrams which contrast the uncloaked scattered pressure amplitude of the spherical shell 
with the scattered pressure amplitude when the cloak is applied. This is done while decreasing 
the transformed radius length to <5 = 0.95a, 5 = 0.5a and 5 = 0.1a. Furthermore, we compute 
the scattering coefficients for these cases and for that of a homogeneous sphere. We were able 
to verify that the acoustic cloak effectively reduces the amplitudes of the pressure forward scat¬ 
tered field. However, as the frequency increases the cloak’s effectiveness tends to decrease, and 
there is often a trade-off between the forward and backscattered amplitudes. An investigation 
into the scattered pressure amplitude in terms of the normal modes revealed greater scattered 
amplitudes by a cloaked shell with a transformed radius length of 5 = 0.95a than by an un¬ 
cloaked sphere, even at lower incident frequencies. However, as shown earlier, this result does 
not indicate net scattered pressure field since the scattered pressure of the various modes are 
not necessarily in phase. Increase frequency results in a shift to higher order modes dominating 
the scattering. This effect is more pronounced when 5 = 0.95, but occurs also when 5 = 0.5 
at 5000 Hz. These results are not completely understood and require further research, including 
an exhaustive investigation throughout the frequency spectrum. 

We also analyzed the backscattered pressure field amplitudes, in terms of the dimensionless 
frequency ka. The effectiveness of the acoustic cloak regarding the backscattered energy has 
particular interest and importance, since most sonar systems are mono-static. We were able 
to verify occurrences of greater scattered amplitude at 5 = 0.95, as well as a shift to higher 
excited mode numbers at higher frequencies. Furthermore, for the backscattered energy, we 
observed that this shifting of the excited modes obeys a consistent relationship proportional to 
the decrease in the transformed radius length. On the other hand, the decrease in 5 corresponded 
to a significant reduction in the backscattered amplitude. Only a reduction to 5 = 0.1 leads to a 
very effective effective cloak throughout all the dimensionless frequency spectrum considered. 
Nonetheless, as the transformed radius is made smaller, acoustic cloak parameters may require 
greater density or lower stiffness values than may be feasible. 

A final analysis considered the acoustic cloak parameters and their behavior within the metafluid. 
We were able to verify that the analytical conditions, imposed with the purpose of making the 
outer surface of the cloak match the sea water physical properties, are satisfied, as we saw all 
the normalized parameters equal unity at r = b, on Figure 4.11. Moreover, we were able to 
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verify that the density parameter in the IC case, and the stiffness parameter in the PM case, 
have an asymptotic behavior at the inner boundary of the cloak, approaching infinity or zero, 
respectively. A reduction in the radius length to 5 = 0.1a would require a density 280 times 
greater than the density of sea water for an IC implementation, or a stiffness very close to zero 
(0.004) for a PM implementation, both of which are unreasonable. Nevertheless, the required 
increase in the density is much greater than the required decrease in stiffness, which leads to 
the conclusion that PM materials constitute a more realistic approach in the fabrication of an 
acoustic cloak. Another alternative that should be explored is the PM-IC metafluid. However, 
it is practically very difficult to design an anisotropic material in density and stiffness. Con¬ 
sequently, careful evaluation of the acoustic cloak parameters must be made and ultimately a 
discrete layered approach would likely be adopted. 
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CHAPTER 6: 

MATLAB MODELING AND SIMULATION 


% / / ************************************************************* 
%// File: Ana_Vieira_thesis.m 
%// Name: Ana Margarida Vieira 

%// MSEA Thesis - Transf. Acoustics applied to Scattering 

%// from a Thin Spherical Shell 

%// Co-Advisor: Professor Brett Borden 

%// Co-Advisor: Professor Clyde Scandrett 

%// Second Reader: Professor Steve Baker 

%// Chair, Engineering Acoustics Academic Committee: 

%// Professor Daphne Kapolka 
%// June 2011 

% / / ************************************************************** 

o, o. 
o o 

clear all; clc; format long 
%% Initialize Variables 
f = 150; % frequency (Hz) 
nth_mode =31; % number of normal modes 
if nth_mode == 0; 

n = 0:nth__mode; % to make possible to see n=0 
nmodes = length(n) ; 

else 

n = 0:nth_mode; 
nmodes = length (n) 

end 

ntheta = 101; 

theta = linspace (0,pi,ntheta); % polar angle theta (0:pi) 
ttheta = linspace (0,2*pi,ntheta); % polar angle theta (0:pi) 
ttheta_deg = ttheta .*180./pi; 

eta = cos(ttheta);% argument of the Legendre Polynomial cos (theta) 

cf = 1500; % sound velocity 

a = 1; % radious of the sphere 

b = 2*a; % outer surface of the cloak 

w = 2*pi*f; % angular frequency (rad/s) 

k = w/cf; % wave number 

f_vec = 0:3000; % range of frequencies 

w_vec = 2*pi*f_vec; % angular frequency (rad/s) 

k_vec = w__vec./cf; 

ka_vec = k_vec.*a; 
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rhof = 1026; % sea water density (kg/m"3) 

bulkf = 2.28*10"9; % sea water isothermal Bulk Modulus (Pa) 
rhos = 7700; % shell density (steel) (kg/m"3) 
nu = 0.28; % Poisson's ratio for the shell 
E = 19.5*10^10; % Young's modulus for the shell (Pa) 

G = 8.3*10 / '10; % shear modulus of the shell 
c_s = sqrt(G/rhos); %shear speed 

lambda_s = c_s/f; % shear wavelength of the shell 
h = .05; % thickness of the shell (m) 

P_i =1; % incident pressure amplitude (m) 
ri = 0.001; 
r = a:ri:l*b; 

kr = k*r; % wave number * distance from the surface of the sphere 

ka = k*a; % wave number*radius at the surface of the sphere 

lambdan = n.*(n+l); % need for the shell radial displacement... 

...amplitude calculations 

cp = sqrt(E/((1-nu"2)*rhos)); 

beta = sqrt(h"2/(12*a"2) ) ; 

alpha = nu + lambdan -1; 

gamma = 1+beta; 

epsilon = 1+nu; 

Omega = (w*a/cp); 
delta = 0.5 
deltal = 0.1; 
delta2 = 0.5; 
delta3 = 0.95; 
hh = num2str(h); 

%% Legendre Polynomials 

for u = 1:length(n); % order nth of the Legendre Polynomial 

Pnm = legendre (u-1, eta); % generates the matrix: rows==... 
...each m from m=0 until m=n; columns==cos(theta) 

PnO (u, :) = Pnm(l, :); % Pick the first row of the matrix p... 

...(so that we get a matrix where each row is ... 

...a vector corresponding to the Legendre Polynomial... 

...of order n and m=0 

PL = PnO; % final matrix where each row is ... 

...the vector corresponding to the Legendre Polynomial... 

...of order n and m=0 

end 

%% Specific Acoustic Impedances of the shell and the fluid 
Zm = (( ( — 11)*rhos*cp)/Omega) * (h/a) *... 

((Omega"4-Omega"2*(gamma*alpha+2*epsilon+beta"2* 
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lambdan . * alpha) - lambdan . * (beta / '2*alpha+epsilon) . "2 + 2* 
epsilon*gamma*alpha+beta"2*alpha."2.*lambdan*gamma))./ 

((Omega"2)-gamma*alpha) ; 

% specific acoustic impedance of the shell 
Za = cf*rhof; % specific acoustic impedance of the fluid 
%% Scatered Pressure Field 
% Near Field 
for v = 1:length(n) 

RnR(v) = -(1i." (v-1) ) .*( (2 *(v-1)+1)*. . . 

( ( ( ( (v-1)/ka)*(besselj(((v-1)+1/2) ,ka)))- 
besselj(((v-1)+3/2),ka)) /((((v-1)/ka)*(besselh(((v-1)+ 

1/2), ka)))-besselh(((v-1)+3/2) , ka) ) ) ) ; 

% nf coefficient rigid be 

RnF(v) = -(li."(v-1)).*(2*(v-1)+1)*((besselj(((v-1)+1/2),ka)... 

/ (besselh(((v-1)+1/2),ka)))) ; 

% nf coefficient "press, release" be 
Rn(v) = -(li."(v-1)).*(2*(v-1)+1)*... 

((Zm(v) *((((v-1)/ka)*(besselj(((v-1)+1/2), ka))) - 
besselj(((v-1)+3/2), ka) ) -... 

(li) *Za *besselj(((v-1)+1/2) , ka) ) / . . . 

(Zm(v) *((((v-1)/ka)*(besselh(((v-1)+1/2), ka)))- 
besselh(((v-1)+3/2),ka))-... 

(li) *Za *besselh(((v-1)+1/2),ka))) ; 

% nf coefficient thin shell without cloak 
Rn_cloak(v) = -(11." (v-1)) .*(2 *(v-1)+1)*... 

((Zm(v) *((((v-1)/k*delta)*(besselj(((v-1)+1/2), 
k*delta)))-besselj(((v-1)+3/2),k*delta))-... 

(li) *Za *besselj(((v-1)+1/2),k*delta) )/. . . 

(Zm(v) *((((v-1)/ka)*(besselh(((v-1)+1/2), 
k*delta)))-besselh(((v-1)+3/2),k*delta))-... 

(li) *Za *besselh(((v-1)+1/2),k*delta))); 

% nf coefficient thin shell without cloak 
% variation of the transformed radius (delta) 

Rn_cloak_deltal(v) = -(li."(v-1)).*(2*(v-1)+1)*... 

((Zm(v) *((((v-1)/k*deltal)*(besselj(((v-1)+1/2), 
k*deltal)))-besselj(((v-1)+3/2),k*deltal) ) + . .. 

(li) *Za *besselj(((v-1)+1/2),k*deltal))/... 

(Zm(v) *((((v-1)/ka)*(besselh(((v-1)+1/2), 
k*deltal)))-besselh(((v-1)+3/2),k*deltal))+... 

(li) *Za *besselh(((v-1)+1/2),k*deltal))) ; 

% nf coefficient thin shell without cloak 
Rn_cloak_delta2(v) = -(li.~(v-1)).*(2*(v-1)+1)*... 
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((Zm(v) *((((v-1)/k*delta2)*(besselj(((v-1)+1/2) , 

k*delta2)))-besselj(((v-1)+3/2),k*delta2) ) + . .. 
(li) *Za *besselj(((v-1)+1/2),k*delta2))/... 

(Zm(v) *((((v-1)/ka)*(besselh(((v-1)+1/2), 

k*delta2)))-besselh(((v-1)+3/2),k*delta2) ) + . .. 
(li) *Za *besselh(((v-1)+1/2),k*delta2))) ; 

% nf coefficient thin shell without cloak 
Rn_cloak_delta3(v) = -(li." (v-1)) .* (2*(v-1)+1)*... 

((Zm(v) *((((v-1)/k*delta3)*(besselj(((v-1) +1/2) , 
k*delta3)))-besselj(((v-1)+3/2),k*delta3) ) + . .. 
(li) *Za *besselj(((v-1)+1/2),k*delta3) )/. . . 

(Zm(v) *((((v-1)/ka)*(besselh(((v-1)+1/2) , 
k*delta3)))-besselh(((v-1)+3/2),k*delta3) ) + . .. 
(li) *Za *besselh(((v-1)+1/2),k*delta3))); 

% nf coefficient thin shell without cloak 
% Far Field (ff) 

RnR_ff(v) = RnR(v)/(li"(v-1)); 

% ff coefficient rigid be 
RnF_ff(v) = RnF(v)/(li"(v-1)); 

% ff coefficient "press, release" be 
Rn_ff(v) = Rn(v)/(li" (v-1) ) ; 

% ff coefficient Shell without cloak 
Rn_cloak_ff(v) = Rn_cloak(v)/ (li"(v-l)); 

% ff coefficient Shell with cloak 

Rn_cloak_deltal_ff(v) = Rn_cloak_deltal(v)/ (li"(v-l)); 

% ff coefficient Shell with cloak 

Rn_cloak_delta2_ff(v) = Rn_cloak_delta2(v)/ (li"(v-l)); 

% ff coefficient Shell with cloak 

Rn_cloak_delta3_ff(v) = Rn_cloak_delta3(v)/ (li"(v-l)); 

% ff coefficient Shell with cloak 

end 

%% Form Function 

amp_polar_RnR_nf = (2/ka).*abs(RnR*PL); 

% form function nf rigid be 
amp_polar_RnF_nf = (2/ka).*abs(RnF*PL); 

% form function nf "press, release" be 
amp_polar_Rn_nf = (2/ka).*abs(Rn*PL); 

% form function nf thin shell without cloak 
amp_polar_Rn_cloak_nf = (2/ka).*abs(Rn_cloak*PL); 

% form function nf thin shell with cloak 
amp_polar_RnR_ff = (2/ka).*abs(RnR_ff*PL); 

% form function ff rigid be 
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amp_polar_RnF_ff = (2/ka).*abs(RnF_ff*PL); 

% form function ff "press, release" be 
amp_polar_Rn_ff = (2/ka).*abs(Rn_ff*PL); 

% form function ff thin shell without cloak 
amp_polar_Rn_cloak_ff = (2/ka).*abs(Rn_cloak_ff*PL); 

% form function ff thin shell with cloak 

amp_polar_Rn_cloak_deltal_ff = (2/ka).*abs(Rn_cloak_deltal_ff*PL); 

% form function ff thin shell with cloak 

amp_polar_Rn_cloak_delta2_ff = (2/ka).*abs(Rn_cloak_delta2_ff*PL); 

% form function ff thin shell with cloak 

amp_polar_Rn_cloak_delta3_ff = (2/ka).*abs(Rn_cloak_delta3_ff*PL); 

% form function ff thin shell with cloak 

%% Scattering Coefficients 

sigcR = -(4/ka~2)*imag(li*sum(RnR_ff)) 

% scattering coefficient rigid be 
sigcF = -(4/ka~2)*imag(li*sum(RnF_ff)) 

% scattering coefficient "press, release" be 
sig_wc = - (4/ka / '2) *imag (li*sum (Rn_ff) ) 

% scattering coefficient spherical thin shell without cloak 
sig_c_deltal = -(4/ka~2)*imag(li*sum(Rn_cloak_deltal_ff)) 

% scattering coefficient spherical thin shell with cloak 
sig_c_delta2 = -(4/ka~2)*imag(li*sum(Rn_cloak_delta2_ff)) 

% scattering coefficient spherical thin shell with cloak 
sig_c_delta3 = -(4/ka~2)*imag(li*sum(Rn_cloak_delta3_ff)) 

% scattering coefficient spherical thin shell with cloak 
%% Backscattering Amplitudes 
for i = 1: length(f_vec) 
for v = 1:length(n) 
k = ka_vec(i); 
ka = k*l; 

RnR(v) = -(li.~(v-1)).*((2*(v-1)+1)*... 

(((((v-1)/ka)*(besselj(((v-1) +1/2) f ka) ) ) - 
besselj ( ( (v-1)+3/2) , ka)) ... 

/ ( ( ( (v-1)/ka)*(besselh(((v-1)+1/2),ka)))- 
besselh ( ((v-1)+3/2) , ka) ) ) ) ; 

% nf coefficient rigid be 

RnF(v) = -(li.~(v-1)).*(2*(v-1)+1)*((besselj(((v-1)+1/2) , ka) 
/ (besselh(((v-1)+1/2) , ka) ) ) ) ; 

% nf coefficient "press, release" be 
Rn_steel(v) = -(li.~(v-1)).*(2*(v-1)+1)*... 

((Zm(v) *((((v-1)/ka)*(besselj(((v-1)+1/2),ka)))- 
besselj(((v-1)+3/2),ka))-. . . 
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(li) *Za *besselj(((v-1)+1/2),ka))/. . . 

(Zm(v) * ( ( ((v-1)/ka)*(besselh(((v-1)+1/2), ka)))- 
besselh(((v-1)+3/2),ka))- . . . 

(li) *Za *besselh(((v-1)+1/2),ka))); 

% nf coefficient thin shell without cloak 
Rn_cloak_deltal(v) = -(li."(v-1)).*(2*(v-1)+1)*... 

((Zm(v) *((((v-1)/k*deltal)*(besselj(((v-1)+1/2), 
k*deltal)))-besselj(((v-1)+3/2),k*deltal))-... 
(li) *Za *besselj(((v-1)+1/2),k*deltal))/... 

(Zm(v) *((((v-1)/ka)*(besselh(((v-1)+1/2) , 
k*deltal)))-besselh(((v-1)+3/2),k*deltal))-... 
(li) *Za *besselh(((v-1)+1/2),k*deltal))); 

% nf coefficient thin shell with cloak (\deltal) 
Rn_cloak_delta2(v) = -(li."(v-1)).*(2*(v-1)+1)*... 

((Zm(v) *((((v-1)/k*delta2)*(besselj(((v-1)+1/2), 
k*delta2)))-besselj(((v-1)+3/2),k*delta2))-... 
(li) *Za *besselj(((v-1)+1/2),k*delta2))/... 

(Zm (v) *((((v-1)/ka)*(besselh(((v-1)+1/2), 
k*delta2)))-besselh(((v-1)+3/2),k*delta2))-... 
(li) *Za *besselh ( ( (v-1)+1/2),k*delta2))); 

% nf coefficient thin shell with cloak (\delta2) 
Rn_cloak_delta3(v) = -(li."(v-1)).*(2*(v-1)+1)*... 

((Zm(v) *((((v-1)/k*delta3)*(besselj(((v-1)+1/2) , 
k*delta3)))-besselj(((v-1)+3/2),k*delta3))-... 
(li) *Za *besselj(((v-1)+1/2),k*delta3))/... 

(Zm(v) *((((v-1)/ka)*(besselh(((v-1)+1/2), 
k*delta3)))-besselh(((v-1)+3/2),k*delta3))-... 
(li) *Za *besselh(((v-1)+1/2) , k*delta3) ) ) ; 

% nf coefficient thin shell with cloak (\delta3) 

% Far Field (ff) 

RnR_f f(v) = RnR(v) /(li" (v-1) ) ; 

% ff coefficient rigid be 
RnF_ff(v) = RnF(v)/(li"(v-1)); 

% ff coefficient "press, release" be 
Rn_steel_ff(v) = Rn_steel(v)/(li"(v-1)); 

% ff coefficient Shell without cloak 

Rn_cloak_deltal_ff(v) = Rn_cloak_deltal(v)/ (li"(v-l)); 

% ff coefficient Shell with cloak 

Rn_cloak_delta2_ff(v) = Rn_cloak_delta2(v)/ (li"(v-l)); 

% ff coefficient Shell with cloak 

Rn_cloak_delta3_ff(v) = Rn_cloak_delta3(v)/ (li"(v-l)); 

% ff coefficient Shell with cloak 
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end 

sumRnR=0; 
sumRnF=0; 
sumRn_steei=0; 
sumRn_cloak_deltai=0; 
sumRn_cloak_delta2=0; 
sumRn_cloak_delta3=0; 
for ii=l:nmodes 

sumRnR=sumRnR+(-1)"(ii—1)*RnR_ff(ii); 

% summing rigid sphere coeffs 

end 

sumRnR=abs(sumRnR); 

amp_polar_RnR_ff_bks(i) = (2/ka_vec(i))*sumRnR; 
for ii=l:nmodes 

sumRnF=sumRnF+(-1)"(ii-1)*RnF_ff(ii); 

% summing "pres, release" sphere coeffs 

end 

sumRnF=abs(sumRnF); 

amp_polar_RnF_ff_bks(i) = (2/ka_vec(i))*sumRnF; 
for ii=l:nmodes 

sumRn_steel=sumRn_steel+(-1)~(ii-1)*Rn_steel_ff(ii); 

% summing shell coeffs (without cloak) 

end 

sumRn_steel=abs(sumRn_steel); 

amp_polar_Rn_steel_ff_bks(i) = (2/ka_vec(i))*sumRn_steel; 
for ii=l:nmodes 

sumRn_cloak_deltal=sumRn_cloak_deltal+(-1)"(ii-1)* 
Rn_cloak_deltal_ff(ii); 

end 

sumRn_cloak_deltal=abs(sumRn_cloak_deltal); 
amp_polar_Rn_cloak_deltal_ff_bks(i) = (2/ka_vec(i))* 
sumRn_cloak_deltal; 
for ii=l:nmodes 

sumRn_cloak_delta2=sumRn_cloak_delta2+(-1)"(ii-1)* 
Rn_cloak_delta2_ff(ii); 

end 

sumRn_cloak_delta2=abs(sumRn_cloak_delta2); 
amp_polar_Rn_cloak_delta2_ff_bks(i) = (2/ka_vec(i))* 
sumRn_cloak_delta2; 
for ii=l:nmodes 

sumRn_cloak_delta3=sumRn_cloak_delta3+(-1)"(ii-1)* 
Rn_cloak_delta3_ff(ii) ; 
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end 

sumRn_cloak_delta3=abs(sumRn_cloak_delta3); 
amp_polar_Rn_cloak_delta3_ff_bks(i) = (2/ka_vec(i))* 
sumRn_cloak_delta3; 

end 

%% Transformation function f(r) 

f_r = delta + (b +a -2 *delta)*((r-a)./(b-a)) + 

(delta - a)*( (r-a) ./(b-a)) ."2; 
f_r_deriv = ( (b+a - 2*delta)/(b-a)) + 

2*(delta -a)*((r-a)/(b-a)) *(1/(b-a)); 

%% Inertial cloak 

rho_r = (f_r_deriv.*((r./f_r).~2))*rhof; 
rho_r_norm = rho_r / rhof; 
rho_r_norm_max = max(rho_r_norm); 
rho_perp = (1./f_r_deriv)*rhof; 
rho_perp_norm = rho_perp / rhof; 
rho_perp_norm_max = max(rho_perp_norm); 
bulk_iso = (1./f_r_deriv).* ((r./f_r).~2)*bulkf; 

bulk_iso_norm = bulk_iso / bulkf; 
bulk_iso_norm_max = max(bulk_iso_norm); 

%% Pentamode material 

rho_iso = f_r_deriv.*(f_r./ r)."2*rhof; 
rho_iso_norm = rho_iso / rhof; 
rho_iso_norm_max = max(rho_iso_norm); 
bulk_r = (1./f_r_deriv).* ((f_r./r)."2)*bulkf; 

bulk_r_norm = bulk_r / bulkf; 
bulk_r_norm_max = max(bulk_r_norm); 
bulk_perp = f_r_deriv.*bulkf; 
bulk_perp_norm = bulk_perp / bulkf; 
bulk_perp_norm_max = max(bulk_perp_norm); 

%% Figures 

% Polar plot rigid BC, "pressure Release" BC, 

% shell without cloak, function of theta 
figure (10) 
if nth__mode ==0 

r_maxl0 = .5; 

else 

r_maxl0 = max(amp_polar_Rn_ff); 

end 

h_fakel0 = polar (ttheta, r_maxlO*ones (size (ttheta) ) , ''w'); 

hold on; 

set(h_fakel0, 'Visible', 'Off'); 
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polar (ttheta, amp_polar_RnR_ff, 'm—') 

hold on 

polar (ttheta, amp_polar_RnF_ff, ' c—') 

hold on 

polar (ttheta, amp_polar_Rn_ff, 'b') 

hold off 
grid on 

legend (' ' , ''Rigid'' , 'Press. Re lease' , ' Sph . She 11','location', 

'NorthEastOutside' ) 
if nth_mode ==0 

title({['Scattered Pressure Field. Far Field. 

Thickness of the Shell ' hh ' ' ] ; 

['Fundamental Mode n = 0']; 

['Frequency ' int2str(f) ' Hz']}) 
else 

title({['Scattered Pressure Field. Far Field. 

Thickness of the Shell ' hh ' ' ] ; 

['' int2str(nmodes-1) ' Mode(s), summed from n = 0 to 

n =' int2str(nmodes-2)'.'];['Frequency' int2str(f) ' Hz']}) 

end 

saveas(gcf, 'polarnocloak', 'pdf') 

% Polar plot. Shell with no cloak and with cloak with 
% variation of the transformed radius (delta) 
figure (14) 
if nth_mode==0 

r_maxl4 = .5; 

else 

r_maxl4 = max(amp_polar_Rn_ff) + 0.1; 

end 

h_fakel4 = polar(ttheta, r_maxl4*ones(size(ttheta)),'w'); 
set(h_fakel4, 'Visible', 'Off') 
hold on 

hi = polar(ttheta,amp_polar_Rn_cloak_deltal_ff, 'm'); 
hold on 

h2 = polar(ttheta,amp_polar_Rn_cloak_delta2_ff, 'r'); 
hold on 

h3 = polar(ttheta,amp_polar_Rn_cloak_delta3_ff, 'g'); 
hold on 

h4 = polar(ttheta,amp_polar_Rn_ff, 'b'); 
hold off 

legend({[' '];['\delta = ' num2str(deltal) '']; 

['\delta =' num2str(delta2) '']; 


65 



r r 


['\delta = ' num2str(delta3) 

['no cloak']},'location','NorthEastOutside') 
if nth__mode ==0 

title({['Scattered Pressure Field. 

Thickness of the Shell ' hh ' ' ] ; 

['Fundamental Mode n = 0'];['Frequency ' int2str(f) ' Hz']}) 

else 

title({['Scattered Pressure Field. 

Thickness of the Shell ' num2str(h) '']; 

['' int2str(nmodes-1) ' Mode(s), summed from n = 0 to 

n =' int2str(nmodes-2)'.'];['Frequency ' int2str(f) ' Hz']}) 

end 

saveas(qcf, 'cloakshell', 'pdf') 

% IC and PM parameters 
figure (30) 

semilogy(r, rho_r_norm, 'm','LineWidth',2) 
hold on 

semilogy(r, rho_perp_norm, 'g','LineWidth',2) 
hold on 

semilogy(r, bulk_iso_norm, 'k-LineWidth',2) 
hold on 

semilogy(r, rho_iso_norm, 'c-LineWidth',2) 
hold on 

semilogy(r, bulk_r_norm, 'r','LineWidth',2) 
hold on 

semilogy (r, bulk_perp_norm, 'b','LineWidth',2) 
grid on 

legend('\rho_r IC','\rho_\perp IC','K IC', 

'\rho PM','K_r PM','K_\perp PM') 
xlabel('Radius of the cloak (a< r <b) ' ) 
ylabel({['Normalized Material Parameters']}) 
if nth_mode==0 

title({['Normalized Material Parameters.']; 

['Thickness of the Shell ' hh '']; 

['Transformed radius (\delta):' num2str(delta) 'a']}) 

else 

title({['Normalized Material Parameters.']; 

['Thickness of the Shell ' hh '']; 

['Transformed radius (\delta):' num2str(delta) 'a']}) 

end 

hold off 

saveas(gcf, 'cart_cloak_params', 'pdf') 
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LineWidth',2) 


% Backscattering energy plot 
figure (40) 

plot(ka_vec,amp_polar_RnR_ff_bks,'m' , ' 
hold on 

plot(ka_vec,amp_polar_RnF_ff_bks,'c' , 'LineWidth',2) 
grid on 

plot(ka_vec, amp_polar_Rn_steel_ff_bks,'b', 'LineWidth',2) 

hold on 

plot(ka_vec,amp_polar_Rn_cloak_deltal_ff_bks,'y','LineWidth',2) 
hold on 

plot(ka_vec,amp_polar_Rn_cloak_delta2_ff_bks,'g','LineWidth',2) 
hold on 

plot(ka_vec,amp_polar_Rn_cloak_delta3_ff_bks,'r','LineWidth',2) 
hold off 
grid on 

legend({['Rigid'];['Press.Release'];['Sph. Shell steel']; 
['Cloak: \delta = ' num2str(deltal) '']; 

['Cloak: \delta = ' num2str(delta2) '']; 

['Cloak: \delta = ' num2str(delta3)];}); 
xlabel('ka') 

ylabel('backscattered pressure amplitude') 

title('Comparison between backscattered amplitudes') 

o. o. 
o o 

%//end of file Ana_Vieira_thesis.m 
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